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Abstract 

A  combined  experimental-numerical  investigation  is  being  conducted  to  study 
dynamic  load  transfer  in  particulate  materials  due  to  explosive  loadings.  The  primary 
emphasis  in  the  study  is  to  relate  the  microstructural  features  of  the  particulate  media  to 
the  load  transfer  process.  The  experimental  technique  of  dynamic  photoelasticity  is  used 
to  investigate  the  effect  of  cementation  and  of  the  particle  shape  on  the  local  contact 
stress  fields.  The  stiffness  of  the  cement  relative  to  that  of  the  particle  controls  the 
location  of  the  peak  contact  stresses.  Strong  cementation  increases  the  load  transfer 
velocity  and  also  promotes  fracture  of  the  particles.  The  particle  shapes  as  presently 
studied  in  this  research  program  seem  to  have  little  influence  on  the  load  transfer  process. 
The  applicability  of  the  fiber  optic  sensors  and  the  speckle  techniques  to  contact  stress 
measurements  is  evaluated.  Fiber  optic  sensors  show  promise  of  future  applications  to 
three  dimensional  problems.  Discrete  element  numerical  wave  simulation  has  been 
conducted  for  saturated  granular  materials  through  the  introduction  of  a  new  contact  law 
using  elastohydrodynamic  theory.  Pore  fluid  acts  to  decrease  the  wave  speed  and  increase 
the  attenuation.  Future  numerical  studies  will  focus  on  additional  changes  of  the 
interparticle  contact  response  through  cementation  and  particle  shape  effects. 
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Statement  of  Work 


The  primary  objective  of  this  investigation  is  to  study  the  role  of  microstructure 
on  the  dynamic  behaviour  of  a  variety  of  particulate  materials.  In  particular,  attention  is 
directed  towards  relating  microstructural  features  including  packing  fabric,  contact 
behaviors,  and  particle  shape  and  size,  to  wave  propagation  characteristics  of  wave  speed, 
amplitude  attenuation,  and  dispersion. 

Status  of  the  Research  Effort 

In  the  past  year  we  initiated  several  different  studies  and  these  include:  (1) 
Relating  granular  fabric  to  wave  propagational  behaviours,  (2)  Evaluating  fiber  optic 
sensor  and  applying  it  to  contact  strain  measurement,  (3)  Applying  white  light  speckle 
photography  to  displacement  measurements,  (4)  Studying  the  effect  of  cementation  on  the 
load  transfer  phenomenon,  (6)  Ljooking  at  the  influence  of  particle  shape  on  the  wave 
propagation  process.  The  details  of  these  efforts  are  presented  in  the  following  chapters. 
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CHAPTER  1 

THE  RELATIONSHIP  OF  GRANULAR  FABRIC  WITH  WAVE 
PROPAGATIONAL  BEHAVIORS 


1.1  Introduction 

Because  wave  propagation  in  granular  materials  is  locally  transmitted  through  and 
across  contact  regions  between  adjacent  particles,  material  microstructure  or  fabric  will  be 
related  to  the  propagational  characteristics  of  wave  speed  and  amplitude  attenuation. 
Fabric  of  granular  materials  can  be  described  as  the  spatial  arrangement  of  the  solid 
particles  and  associated  voids.  Over  the  past  several  years,  many  fabric  measures  have 
been  proposed  for  granular  media,  but  no  unique  choice  of  fabric  description  has  been 
universally  accepted  which  can  adequately  describe  the  mechanical  response  of  such 
media. 

Some  research  has  categorized  fabric  into  two  types:  orientation  fabric 
(orientation  of  individual  particles)  and  packing  fabric  (mutual  relation  of  individual 
particles  to  the  others).  Orientation  fabric  may  be  quantitatively  defined  by  a  vector  mean 
direction  and  a  vector  magnitude  to  characterize  orientation  of  non-spherical  or  non¬ 
circular  particles.  This  fabric  measure  is  commonally  represented  by  an  angular  measure 
(with  respect  to  a  reference  direction)  of  the  long  axes  of  individual  particles.  Packing 
fabric  measures  have  included  branch  vectors,  normal  contact  vectors,  coordination  or 
contact  numbers,  void  characteristics,  etc.  Examples  of  these  fabric  measures  are 
illustrated  schematically  in  Fig.  1.1.  Some  of  these  fabric  measures  are  kinematical  in 
nature,  determined  primarily  by  the  particle  shape  and  packing  geometry.  Other  measures 
such  as  those  related  to  the  contact  conditions  are  kinetic,  and  are  detenu  ined  by  the 
particle  material  and  shape  properties  and  also  by  the  contact  surface  condii.ons. 

Fig.  1.2  shows  a  photograph  of  the  dynamic  photoelastic  fringe  patterns  associated 
with  an  actual  wave  moving  through  a  model  particulate  medium  containing  scmie  of  the.se 
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fabric  structures.  It  has  been  observed  for  both  static  and  dynamic  loading  conditions  that 
particulate  materials  transmit  mechanical  loadings  along  a  series  of  complex  discrete  paths. 
These  discrete  paths  are  established,  in  relation  to  the  loading  direction,  by  many  of  the 
microstructural  variables  mentioned,  and  therefore  local  wave  propagation  is  determined 
by  granular  fabric  through  the  creation  of  local  wave  guides. 

This  section  outlines  our  research  which  examines  the  relationship  of  granular 
fabric  with  the  wave  propagational  variables.  The  model  granular  materials  under  study 
will  be  large  random  assemblies  of  circular  particles  which  have  been  created  using  one  or 
more  of  our  random  media  generator  codes.  This  study  focussed  primarily  on  the  effects 
of  packing  fabric,  including  branch  vector  distributions  (which  for  circular  particles 
coincide  with  normal  contact  vectors),  path  fabric,  and  void  polygon  fabric. 

1.2  Branch  Vector  Fabric 

Branch  vectors  have  been  used  as  fabric  measures  in  several  previous  static 
investigations  of  granular  media,  and  this  measure  also  looks  appropriate  for  wave 
propagation  problems  since  dynamic  loads  are  transmitted  primarily  through  contact 
points  which  lie  along  branch  vectors.  Thus  it  seems  reasonable  to  assume  that  branch 
vector  distributions  could  be  related  to  the  wave  propagational  behaviors  in  such  media. 

To  investigate  such  relationship.s,  .several  of  the  large  random  assembly  generators 
previosly  discussed  in  an  earlier  AFOSR  report  (Shukla  and  Sadd,  1992)  were  used  to 
construct  six  assemblies  ’.vi;h  dit'terenl  microstructures.  These  assemblies  are  named  as  S- 
1,  M'l.  M-2,  W-1,  W-2.  and  W-3.  Shown  in  Fig.  1.3,  assembly  S-1  was  generated  by  the 
strongly  anisotropic  generator.  Assemblies  M-1  and  M-2,  displayed  in  Figs.  1.4  and  1.5. 
were  constructed  by  the  moderately  anisotropic  generator.  Finally,  the  weakly  anisotropic 
generaU)r  was  used  to  constructe  as.semblies  W-1,  W-2  ,  and  W-3,  and  they  are  slunvn  in 
Figs  1.6  -  1.8.  The  total  numbers  of  particles  for  each  of  the  six  assemblies  are  listed  in 
Table  1.1.  The  void  ratio  (volume  of  void  /  volume  of  particle)  and  average  coordination 
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number  N  (average  contact  number  per  particle)  are  also  calculated  and  are  given  in  Table 
1.1  for  the  six  assemblies.  Assemblies  S-1  has  the  highest  void  ratio  of  0.43  and  lowest 
coordinate  number  of  2.87,  while  W-3  has  the  lowest  void  ratio  of  0.22  and  the  highest 
coordinate  number  of  4.30. 

The  branch  vectors  of  these  assemblies  were  then  calculated  and  normalized  to 
construct  unit  vectors.  To  find  a  load  preferred  direction,  the  x-component  fc^and  y- 
component  of  each  branch  vector  b  were  determined  and  summed  over  the  entire 
assembly.  Denoting  F*"  as  the  ratio  of  the  summed  y-components  over  the  x-components 
gives 

(1-1) 

The  values  of  F*  for  the  six  assemblies  are  listed  in  Table  1.1.  The  branch  vector  angular 
distributions  (rose  diagrams)  of  the  six  assemblies  are  shown  in  Figs.  1.9  -  1.14.  For 
assembly  S-1  (see  Fig.  1.9),  all  branch  vectors  are  concentrated  in  the  region  defined 
by  approximately  ±30°  from  the  vertical  direction.  No  branch  vector  is  found  along  the 
horizontal  direction,  and  is  2.06  for  this  assembly.  The  branch  vector  distributions  for 
M-1  and  M-2  are  shown  in  Fig.  1.10  and  Fig.  1.11,  respectively.  The  branch  vectors  lie 
mainly  in  vertical  and  horizontal  directions  with  preference  in  vertical  direction,  and  both 
assemblies  have  the  same  F^  ratio  of  1.19.  Figs.  1.12  -  1.14  show  the  branch  vector 
distributions  for  as.semblies  W-1,  W-2,  and  W-3.  The  weakly  anisotropic  assemblies  do 
n 'I  shown  any  significant  preferred  branch  vector  distributions.  The  F'’  ratios  for  these 
cases  are  1.01,  0.9,  and  1.01,  re.spectively.  Using  the  discrete  element  method  with  a 
mmlinear  hysteretic  contact  law,  our  wave  propagation  code  WAPRIPM  (Wave 
P-  -Mgation  in  Particulate  Media)  was  used  to  simulate  the  wave  propagational  behaviors 
(wave  speed  and  amplitude  attenuation)  in  the  six  assemblies  under  study.  The 
siinulatUms  involved  the  comparison  of  the  propagation  of  plane  type  waves  nnwing  along 
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orthogonal  directions  (horizontal  and  vertical)  in  the  generated  assemblies.  The  transient 
input  loading  was  modeled  using  a  triangular  pulse  with  a  peak  value  of  1000  N  and  a 
period  of  60  /is.  To  represent  a  planar  input  wave,  the  time-dependent  pulses  were 
simultaneously  applied  to  particles  along  one  of  the  horizontal  or  vertical  boundaries  of 
the  assemblies  under  study.  To  calculate  the  transmitted  or  output  wave  pulse,  an 
imaginary  horizontal  or  vertical  line  was  drawn  near  the  boundary  opposite  to  where  the 
input  loadings  were  applied.  If  a  branch  vector  of  a  pair  of  particles  in  contact  or  potential 

contact  is  intercepted  by  this  imaginary  line,  the  normal  contact  load  component 
perpendicular  to  the  imaginary  line,  i.e.  either  F,  or  F^,  was  recorded.  In  this  fashion,  the 

contribution  of  the  individual  particle  contacts  could  be  determined,  and  these  recorded 
loads  were  then  summed  and  normalized  with  respect  to  the  sum  of  the  peak  values  of  the 
input  loadings.  These  normalized  contact  loads  are  given  by 


Fyyit) 


Fpeakinpul 

Vf 

^  *  peak.inpul 


(1.2) 


where  the  summation  in  the  numerator  is  over  the  number  of  branch  vectors  intercepted  by 
the  imaginary  line,  while  the  summation  in  the  denominator  is  over  the  number  of  input 
loadings. 

Fig.  1.15  shows  this  wave  propagation  simulation  for  the  media  model  S-1.  The 
normalized  normal  contact  load  component  profiles  F„(/)  and  F„,(/)  corresfnmd  tt) 


vertical  and  horizontal  input  loadings  respectively,  and  the  ratios  of  the  peak  values  of  the 
profiles,  denoted  by 


(1.3) 
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are  given  in  Table  1.1. 

For  assembly  S-1,  is  200.68.  It  is  obvious  that  a  wave  propagating  along  the 

vertical  direction  has  much  less  attenuation  than  one  traveling  along  the  horizontal 
direction.  Also  note  that  the  wave  speed  as  determined  by  the  arrival  time  of  the  averaged 
profile  is  different  for  the  two  propagational  directioizs,  and  the  ratio  of  vertical  wave 
speed  over  horizontal  wave  speed  is  approximately  3.  Note  that  results  for  the  horizontal 
wave  are  not  visible  because  of  the  scale  used  in  Fig.  1.15.  These  results  follow 
qualitatively  from  the  local  branch  vector  distribution  diagram  in  Fig.  1.9.  Since  most 
branch  vectors  lie  along  the  vertical  direction,  the  load  transfer  paths  in  that  direction  are 
rather  continuous  and  straight,  while  just  the  opposite  would  be  true  for  load  transfer 
paths  along  the  horizontal  direction.  Therefore,  assembly  S-1  is  an  anisotropic  medium, 
and  demonstrates  high  directionally  dependent  wave  propagation. 

The  normalized  normal  contact  load  profiles  for  model  assemblies  M-1  and  M-2 
are  shown  in  Figs.  1.16  and  1.17.  For  these  assemblies,  the  discrete  element  simulations 
indicate  that  the  vertical  wave  attenuation  is  still  smaller  than  horizontal  wave,  but  these 
differences  are  not  as  pronounced  as  in  assembly  S-1.  The  wave  traveling  along  vertical 
direction  was  faster  than  the  horizontal  wave,  but  differences  between  the  two  wave 
velocities  is  not  as  large  as  that  for  the  S-I  assembly.  Therefore,  the  degree  of  anisotropy 
is  larger  in  assembly  S-1  than  in  M-1  or  M-2,  which  agrees  with  their  re.speciive  F''  values. 
However,  the  F*  ratio  fails  to  correlate  with  the  anisr>tropic  wave  propagaticmal  behaviors 
of  assemblies  M-1  and  M-2.  The  ratios  of  the  peak  values  of  the  profiles  F  are  33.33 

and  100.02  for  M-1  and  M-2,  re.spectively,  and  the  ratio  of  vertical  v  avc  speed  to 
horizontal  wave  speed  for  M-2  (1.60)  is  greater  than  that  for  M-1  (1.52).  All  of  these 
results  indicate  that  assembly  M-2  has  higher  wave  propagaticmai  amsoiropy,  but 
unfortunately  the  fabric  measure  F’  of  M-2  has  the  same  value  as  M- 1 .  'flu  \  eilical  wave 
profile  in  Fig.  1.17  has  a  small  precursor  pulse  before  the  main  profile.  The  reason  for  this 
is  that  the  anisotropy  within  whole  asserr.bly  is  not  uniform,  and  a  .small  portion  of  the 
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total  wave  will  propagate  along  short  paths  and  will  thus  arrive  sli^tly  ahead  of  the  main 

signal.  The  simulation  results  for  the  weakly  anisotropic  assemblies  are  shown  in  Figs. 
1.18  -  1.20.  The  ratios  for  these  assemblies  are  much  smaller  than  for  assemblies  M-1 

and  M-2.  They  are  1.76  for  W-1, 0.94  for  W-2,  and  1.23  for  W-3.  Compared  with  the  S- 
and  M-assemblies  previously  discussed,  differences  between  vertical  and  horizontal  wave 
speeds  are  relatively  small  for  the  W-assemblies.  This  agrees  with  the  local  branch  vector 
distribution  diagrams  which  show  that  the  W-assemblies  are  weakly  anisotropic.  Again, 
the  ratio  F*’  fails  to  provide  a  correlation  with  fabric  and  wave  behaviors  for  the  W- 
assemblies. 


1.3  Path  Fabric 

When  waves  propagate  through  granular  media,  the  local  or  micro-dynamic  loads 
are  never  uniformly  transmitted  across  a  sample,  but  rather  are  transferred  through  specific 
chains  of  particles  linked  through  contact.  Fig.  1.21  illustrates  the  contact  load 
distribution  (the  thickness  of  the  lines  scaled  to  the  current  maximum  force)  in  a  particle 
bed  due  to  particle  impact  at  the  top  (Thornton  and  Randall,  1988).  Between  these 
special  chains,  particles  may  carry  little  or  no  load.  It  appears  that,  by  chance,  some 
sequences  of  particles  offer  contiguous  paths  through  the  assembly  that  are  straighter  than 
other  possible  paths.  It  is  these  paths  that  carry  most  of  the  load  transfer  rather  than  their 
more  tortuous  neighbors.  Thus,  within  an  assembly,  if  in  one  direction  there  are  more 
contiguous  and  straight  iiaths  than  in  other  directions,  the  wave  should  propagate  more 
easily  in  the  former  direction.  This  section  describes  a  study  of  this  path  fabric  effect  on 
wave  propagation. 

A  path  is  defined  as  a  set  of  continuous  branch  vectors  of  particles  in  contact,  and 
it  was  olxserved  expennientally  by  Shukla  et.al.  (1988)  that  dynamic  load  can  be 
transferred  along  the  path  only  when  every  pair  of  neighlxuing  branch  vectors  b-  and  bj 

in  the  path  satisfy  the  relaiion.ship 
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(1.4) 

which  means  that  the  branch  angle  (p  shown  in  Fig.  1.22  is  less  than  90°.  Therefore,  a 
transferable  or  propagator  load  path  can  be  defined  as  a  series  of  particles  which  are  in 
contact  with  each  other,  such  that  the  dot  product  of  any  pair  of  neighboring  branch 
vectors  is  positive. 

Within  this  section,  two  different  methods  of  constructing  a  path  fabric  measure, 
based  on  the  straightness  of  the  path,  will  be  discussed.  In  one  method,  the  dot  products 
of  all  adjacent  branch  vectors  in  a  given  path  are  calculated,  and  they  are  used  as  basic 
units  to  build  a  path  fabric  measure.  For  example,  a  path  fabric  between  two  arbitrary 
particles  can  be  expressed  as 


(1.5) 


where  inner  summation  is  over  a  path  connecting  the  particles  and  outer  summation  is 
over  all  of  the  possible  paths.  For  the  straight  chain  case,  there  exists  only  one  path 


between  the  two  particles,  and  equation  (1.5)  gives  F’  =  - ,  where  N  is  the  total 

N  -1 


number  of  particles  in  the  path.  On  the  other  hand,  if  there  is  only  a  single  path  between 
two  particles  in  a  two  dimensional  assembly,  and  if  any  branch  angle  in  this  path  equals  90 
degrees,  equation  (1.5)  gives  =  0  .  The  path  will  thus  become  a  non- propagator  if 
/>  •  bj  s  0.  and  in  the  following  sections,  only  propagator  paths  will  be  included  in  the 

investigation. 

A  second  method  to  build  a  path  fabric  measure  uses  x  and  y  branch  vector 
ccmiponents  instead  of  the  dot  product  used  in  equation  (1.5).  One  fabric  measure 
i'ciween  two  points  can  be  written  as 
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..IM 
'  22- 

.  22hk' 
’■  22* 


(1.6) 


where  outer  summations  are  over  all  possible  paths  between  the  two  points,  and  inner 
summations  are  over  individual  paths.  Equations  (1.6)j  and  (1.6)2  are  used  to  describe  the 
path  fabric  along  the  x  and  y  directions,  respectively.  For  granular  media  simulated  by 
circular  disks  of  the  same  size,  the  value  of  equations  (1.6)j  and  (1.6)2  will  not  be  larger 
than  the  disk  diameter,  and  will  equal  to  the  diameter  only  when  the  path  is  a  straight  chain 
along  the  x  or  y  direction. 

Using  the  first  proposed  path  fabric  measure  to  analyze  the  microstructure  of  an 
assembly,  it  is  required  to  find  every  path  between  each  pair  of  boundary  particles,  of 
which  the  input  loading  is  applied  at  one  and  output  contact  load  is  recorded  at  the  other. 
Since  two  paths  are  defined  as  different  if  one  path  includes  any  particle  the  other  does 
not,  the  number  of  the  paths  across  an  assembly  will  be  very  large.  Theoretically,  all 

different  paths  can  be  found  and  ^  ^  can  be  calculated  for  each  path.  However,  for 

large  assemblies  the  total  number  of  paths  between  two  boundary  points  approach  10^  - 
10®,  and  considerable  CPU  time  and  large  memory  space  will  be  needed  to  account  for 
each  path.  For  example,  to  calculate  the  path  fabric  measure  f  between  two  points 
which  lie  on  opposite  boundaries  of  assembly  W-1,  takes  more  than  1  hour  of  CPU  time 
(VAXstation  3100).  Therefore  it  appears  impractical  to  use  equation  1.5). 

One  solution  to  this  problem  is  to  divide  the  whole  assembly  into  several 
subregions  and  use  the  fabric  measure  of  equation  (1.6).  Using  this  scheme  a  general 
assembly 
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can  be  divided  into  subparts,  and  Fig.  1.23  shows  four  such  subregions  for  the  calculation 
of  vertical  paths.  In  this  figure,  dots  stand  for  the  contact  points,  and  the  solid  polygon  line 
represents  one  of  many  paths,  which  start  at  point  A',  go  through  points  B',  C,  and  D',  and 
•  end  at  point  E'.  Equation  (1.6)  for  the  vertical  paths  from  A'  to  E'  can  be  written  as 


Ff 


f; 


(1.7) 


where  |^'^j  is  the  vector  from  A'  to  E'.  The  summation  in  the  numerator  is  over  all 
possible  paths  from  A'  to  E'.  The  equation  can  be  further  rewritten  as 


fP  «  path-A'E- 

..I.  ^ 


p€trticle 


JL 

paih‘A'E' 

p<ut^A'P 


(1.8) 


where  is  the  total  number  of  paths  between  A'  and  E'  ,  and  is  the  total 

number  of  particles  along  a  path.  If  the  paths  are  restricted  through  points  B',  C,  and  D', 
the  total  number  of  paths  is  given  by 


N  =  N  •  N  ■  N  •  N 

path- A- r  path- A' If  path-BC  palh-Cl  r  path- D 


(1.9) 


and  the  total  number  of  particles  equals 


10 


^ parliel»-A'P  “  (  / 

'  ^piTticu)  '  ^patk-ffC  '  ^ patk-Ciy  '  ^ palh-CTP 

A'Bt 

+  (  ' 

^  ^parlicu)  ‘  ^ path- A' ET  '  ^  path-CO  '  ^  path- If  P 

pan 

rtC 

+  (  ' 
A 

^  ^ partial^  '  ^ path- A' Bt  '  ^  path- BC  '  ^ path- DTP 

pail 

rro- 

+  (  ' 
A 

y  ^ptrlicle)  ’  ^ path- A’B  '  ^ path- BC  '  ^ path-CD' 

pail 

In  Fig.  1.23,  taking  the  bottom  subregion  as  an  example,  all  vertical  i  rom  A' 
to  B'  can  be  found,  and  the  particle  numbers  of  all  the  paths  are  summed.  The  above 
procedure  is  repeated  for  the  other  three  regions,  and  equation  (l.S)^  gives  the  path  fabric 
measure  for  all  the  paths  between  A'  and  E’  which  pass  through  B',  C,  and  D'.  Now  let  B' 
take  all  points  on  level  B-B,  C  take  all  points  on  level  C  -  C,  and  D’  take  all  points  on  level 
D-D,  then  equation  (1.8)2  gives  the  fabric  measure  for  all  paths  from  A'  to  E'.  The 
path  fabric  measure  for  the  whole  assembly  can  then  be  obtained  by  summing  the  values 
given  by  equation  (1.8)2  paths  between  boundaries  A-  A  and  E  -  E.  Similarly, 

the  horizontal  value  can  be  calculated  from  equation  (1.8)^.  It  should  be  pointed  out  that 
when  a  vertical  path  is  under  consideration,  the  path  will  be  discarded  if  it  runs  beyond  a 
vertical  zone  of  fixed  width  centered  at  the  initiating  particle  within  a  subregion.  This  is 
reasonable  since  a  highly  zigzagged  path  will  transmit  negligible  load,  and  by  discarding 
these  paths  a  reduction  of  computational  effort  will  result. 

Using  equation  (1.8),  the  vertical  and  horizontal  path  fabrics  are  computed  for  the 
six  assemblies  discus.sed  in  the  previous  section.  The  ratios  of  the  summation  of  vertical 
fabric  to  the  summation  of  hiuizontal  fabric 


= 


vernal 


horizontal 


(1.11) 
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are  given  in  table  1.1.  For  the  strongly  anisotropic  assembly  S-1,  there  is  no  contiguous 
path  along  horizontal  direction  and  this  gives  a  ratio  of  infinity,  which  correlates  with  the 
discrete  element  results  which  show  that  the  horizontal  wave  is  almost  blocked.  The 
ratio  has  the  values  of  5.88  and  2.08  for  the  two  moderately  anisotropic  assemblies  M-1 
and  M-2.  Thus,  there  are  more  paths  along  the  vertical  direction  rather  than  the  horizontal 
direction,  and  waves  prefer  to  propagate  along  the  vertical  direction  in  these  assemblies  as 
predicted  by  discrete  element  analysis.  However,  as  with  the  branch  vector  case,  the  path 
fabric  also  fails  to  indicate  the  stronger  anisotropy  of  assembly  M-2  over  M-1.  The  reason 
for  this  lack  of  correlation  may  be  due  to  gap  closing,  which  can  change  the  path  fabric 
significantly.  Taking  assembly  M-2  as  an  example,  during  the  process  of  vertical  wave 
propagation,  it  appears  that  more  than  40  new  contacts  have  been  created.  The  F’’  ratios 
of  the  weakly  anisotropic  assemblies  W-1,  W-2,  and  W-3  are  0.95,  1.00,  and  1,16, 
respectively,  which  indicates  that  these  assemblies  are  indeed  weakly  anisotropic. 
However,  the  F"  ratio  fails  to  provide  a  correlation  with  path  fabric  and  wave  behaviors 
for  the  W-assemblies. 


1.4  Void  Polygon  Vector  Fabric 

In  the  previous  two  sections,  fabric  measures  associated  with  particles  were 
discussed.  However,  besides  the  solid  particle  phase,  a  granular  material  also  consists  of 
a  void  phase,  and  thus  it  would  seem  reasonable  that  a  fabric  measure  based  on  voids 
could  be  useful  to  correlate  with  wave  propagation. 

To  describe  a  void  with  N  curved  segments,  a  polygon  is  used  as  shown  in  Fig. 
1 .24.  The  polygon  consists  of  the  branch  vectors  linking  the  particles  around  a  void,  and 
thus  a  void  polygon  represents  the  relations  of  this  special  group  branch  vectors. 
A.ssuming  that  a  void  surrounded  by  N  particles  has  N  branch  vectors  bj,  and  b^  a 
tensor  can  be  defined  (similar  to  Konishi  and  Naruse,  1988)  as 
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(1.12) 


Pa  - 

where  and  are  x  and  y  components  of  bj^  respectively.  The  principal  values  of  are 


A. f>*  cos28»  j  +  sin  26, 
The  difference  between  the  principal  values  is  denoted  by 

H  ~  Pi -pi 

and  the  major  principal  direction  is  calculated  as 


Y 

j 


(1.13) 


(1.14) 


0,  -  — arctan 
'  2 


r  N  ^ 

■^sin  20^ 


O'*cos20* 


(1.15) 


vra  j 

Now  the  void  can  be  described  by  Pi,P2,  and  0^:  //  -  pj  -  Pj  is  related  to  the  anisotropy 
of  the  void,  S  =  pjPj  represents  the  void  area,  and  0^  gives  the  void  direction. 

Konishi  and  Naruse  (1988)  constructed  a  local  void  vector  p  whose  value  is  H  and 
direction  is  given  by  0^.  In  this  study,  the  local  void  vector  p  is  defined  as  a  vector 

whose  magnitude  is  HS .  The  reason  for  including  S  in  the  void  vector  definition  is  that 
when  waves  propagate  through  a  granular  material,  the  possibility  of  voids  and  gaps 
closing  is  related  to  void  size.  The  void  vector  direction  will  also  use  the  measure  0^. 

Therefore,  voids  within  a  granular  mass  can  be  found  and  the  void  vectors  may  be 
constructed  for  each  of  them. 

Using  this  scheme  of  local  void  vectors,  a  rose  diagram  similar  to  the  branch 
vector  plots  in  Figs.  1.9  -  1.14  can  be  plotted.  Figs.  1.25  -  1.30  show  void  vector 
distributions  for  the  six  as.semblies  discussed  in  the  previous  sections.  An  examination  of 
these  plots  along  with  the  corre.sponding  branch  vector  diagrams  reveals  that  void  vector 
reflects  the  microstructure  more  clearly  than  branch  vectors.  In  the  void  vector 
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distribution  diagram  for  assembly  S-1,  all  the  void  vectors  lie  in  vertical  direction  to  block 
horizontal  wave  paths,  while  the  corresponding  branch  vector  diagram  fails  to  demonstrate 
this.  For  assemblies  M-1  and  M-2,  void  vector  rose  diagrams  show  an  overwhelming 
majority  of  void  vectors  lying  in  the  vertical  direction.  This  agrees  with  the  numerical 
simulations,  while  the  corresponding  branch  vector  diagrams  only  show  slightly  more 
branch  vectors  in  vertical  direction  than  in  the  horizontal.  The  void  vector  rose  diagrams 
for  the  W-assemblies  also  indicate  preferred  wave  propagation  directions,  which  was  not 
apparent  in  the  branch  vector  diagrams. 

The  summation  of  all  void  vector  x  and  y  components  can  be  given  by  ^  and 
^  respectively,  and  thus  we  can  define  a  ratio 


F'' 


Ip. 


(1.16) 


which  can  be  used  to  measure  the  anisotropy  of  a  medium.  The  values  of  F*”  for  the  sbc 
assemblies  are  given  in  Table  1.1.  This  ratio  increases  with  the  degree  of  assembly 
anisotropy  except  for  assemblies  W-2  and  W-3,  and  it  obtains  a  maximum  value  of  85.08 
in  assembly  S-1,  indicating  a  strongly  anisotropic  medium.  For  moderately  anisotropic 
assemblies  M-1  and  M-2,  the  ratio  decreases  to  2.44  and  2.55  respectively,  which 
successfully  predicts  a  stronger  anisotropy  for  assembly  M-2  over  M-1.  The  values  of  F' 
for  the  weakly  anisotropic  assemblies  W-1,  W-2  and  W-3  are  1.12,  1.05,  and  1.00, 
respectively.  F"  correctly  indicates  that  the  degrees  of  ani.sotropy  in  W-2  and  W-3  are 
weaker  than  in  W-1,  but  as  with  the  other  fabric  measures  it  fails  to  pc^int  luit  W-3  has  the 
smallest  degree  of  anisotropy  among  the  W-assemblies.  Thus,  it  can  be  concluded  that  the 
void  fabric  measure  correlates  with  the  wave  propagational  behaviors  (4'  particulate 
materials  better  than  the  branch  vector  or/and  path  fabric  schemes  discussed  in  the 
previous  sections. 
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1  Assembly 

Particle  Number 

Void  Ratio 

N 

F. 

m 

F" 

n 

B 

822 

0.43 

2.87 

200.68 

2.06 

no  x-path 

85.08 

M-1 

778 

0.28 

3.43 

33.33 

1.19 

5.88 

2.44 

M-2 

854 

0.28 

3.41 

100.02 

1.19 

2.08 

2.58 

W-1 

1042 

0.25 

4.17 

1.76 

1.01 

0.95 

1.12 

W.2 

1296 

0.25 

4.04 

0.94 

0.99 

1.00 

1.05 

W-3 

1338 

0.22 

4.30 

1.23 

1.01 

1.16 

1.00 

Table  1.1  Various  Fabric  Measures  with  Density  and  Coordinate  Number  for  Six  Assemblies 
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Contact 

Conditions 


Particle  Orientation 

Contact  Normal 
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Particle  Density 


Load  Path  Geometry 


and  Elastic  Modulus 


Fig.  1. 1  Common  Fabric  Me^ures  for  Granular  Materials 


I'ig  1 .2  Photoelastic  Fringe  Patterns  of  Wave  Propagation  in 
a  Model  Particulate  Media 
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1-ig  1  4  A  Moderately  Anisotropic  Granular  Assembly  (M- ' ) 
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Fig.  1.5  A  Moderately  Anisotropic  Granular  Assembly  (M-2) 
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Fig.  1.11  Branch  Vector  Distribution  of  Assembly  M-2 
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Fig.  1.17  Normalized  Contact  Load  vs  Time  for  the  Assembly  Shown  in  Fig.  1 .5 
NORMALIZED  CONTACT  LOAD  VS.  TIME 
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Fig.  1.18  Normalized  Contact  Load  vs  Time  for  the  Assembly  Shown  in  F  ig.  1 .6 
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NORMALIZED  CONTACT  LOAD  VS.  TIME 


TIME  (MICROSECOND) 

Fig.  1.19  Nonnalized  Contact  Load  vs  Time  for  the  Assembly  Shown  in  Fig.  1.7 
NORMALIZED  CONTACT  LOAD  VS,  TIME 
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Fig.  1  20  Normalized  Contact  Load  vs  Time  for  the  .Assembly  Shown  in  Fig.  1 .8 
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Fig.  1.25  Void  Vector  Distribution  of  Assembly  S-1 


l  ijj,  !  2('i  Void  Vecloi  Disiribuiion  of  Assembly  M- 1 
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Fig.  1.27  Void  Vector  Distribution  of  Assembly  M-2 


Fig  1  28  Void  Vecloi  Distribulion  <>!  \>,scnihl\ W  i 
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CHAPTER  2 


EVALUATION  OF  A  FIBER  OPTIC  SENSOR  FOR  STRAIN 
MEASUREMENTAND  PRELIMINARY  APPLICATIONS 
TO  CONTACT  MECHANICS 


2.1  Introduction 

In  this  chapter,  an  attempt  has  been  made  to  study  the  applicability  of  both 
attached  and  embedded  fiber  optic  sensors  to  measure  in  plane  and  out  of  plane  strain  due 
to  contact  loadings.  Both  two  and  three  dimensional  experiments  have  been  conducted. 
Prior  to  applying  fiber  optic  sensors  to  contact  mechanics  problems,  they  were  evaluated 
to  determine  the  maximum  strains  that  could  be  measured  reliably.  In  the  past,  resistive 
foil  strain  gages  have  been  used  successfully  to  evaluate  contacts  in  both  static  and 
dynamic  problems  (Xu  and  Shukla,  1993).  However,  strain  gages  are  limited  in  their  use 
to  surface  strain  measurements.  Also,  they  cannot  be  easily  bonded  to  certain  kinds  of 
materials  like  rock  and  concrete.  Fiber  optic  sensors  on  the  other  hand  can  be  used  for 
both  surface  (Sirkis  and  Taylor,  1988)  and  interior  strain  measurements  (Murphy  el.  ai, 
1989).  Fiber  optic  sensors  have  also  been  shown  to  be  applicable  for  use  in  electrically 
noisy  environments  (Griffiths,  1991),  and  in  high  temperature  areas  (Wang,  1992).  The 
ultimate  aim  of  the  study  is  to  accurately  predict  local  strains  at  any  point  in  a  three 
dimensional  body  for  both  static  and  dynamic  loading. 


2.2  Sensor  Construction 

2.2,1  Mach-Zehnder  Interferometer 
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Construction  of  a  Mach-Zehnder  sensor  is  quite  simple,  and  a  schematic  of  the 
setup  is  shown  in  Fig.  2.1.  A  laser  beam  is  split  with  a  beam  splitter  and  the  two  resulting 
beams  are  then  coupled  into  two  separate  single  mode  optical  fibers.  One  of  the  fibers  is 
chosen  as  the  reference  path  and  the  other  fiber  becomes  the  sensing  arm.  The  two 
separate  beams  are  then  recombined  and  the  resulting  interference  pattern  is  projected 
onto  a  photodiode.  The  interference  fringes  produced  can  be  related  to  the  axial  strain 
seen  by  the  sensor. 

In  operation,  the  fringe  pattern  remains  stationary  when  the  sensor  undergoes  no 
strain  and  the  output  from  the  photodiode  is  constant.  When  the  sensor  experiences  axial 
strain,  the  fringe  pattern  shifts  across  the  face  of  the  photodiode  resulting  in  a  sinusoidal 
type  output.  It  is  this  sinusoidal  output  that  can  be  related  to  the  applied  strain. 

Care  must  be  taken  with  the  Mach-Zehnder  interferometer  because  all  of  the  fiber 
in  the  setup  has  the  potential  to  affect  the  light  propagating  through  it.  In  other  words,  the 
entire  length  of  fiber,  in  either  arm,  can  be  the  sensing  element.  Care  must  be  taken  to 
ensure  that  only  the  desired  measurand  is  affecting  the  output  signal.  This  is  one  of  the 
drawbacks  to  the  sensor  and  therefore,  it  is  not  easily  applicable  to  field  applications. 

2.2.2  Fabry-Perot  Interferometer 

Fig.  2.2  shows  a  .schematic  of  an  extrinsic  Fabry-Perot  sen.sor.  1  he  experimental 
setup  is  almo.st  the  .same  as  that  for  the  Mach-Zehnder  and  is  shown  in  Fig.  2.3.  The 
schematic  of  the  actual  sensor.  Fig.  2.2,  shows  two  pieces  of  optical  fihei  tiiat  are  bonded 
into  a  glass  tube.  The  fiber  ends  are  partially  mirrored  to  provide  the  necessary  reflections 
for  a  low  finesse  Fabry-Perot  cavity.  Fiber  diameter  that  is  u.sed  in  the  laboratory  is  125 
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microns  due  to  the  availability  of  the  appropriate  couplers  and  equipment  for  this  size 
fiber.  Purchasing  ready  made  sensors  is  quite  expensive  and  since  the  sensors  may 
become  damaged  during  experiments,  a  simple  method  to  manufacture  them  was  needed. 

Since  this  sensor  is  an  interferometric  sensor,  a  single  mode  fiber  rated  to  carry 
the  particular  input  wavelength  is  used.  A  helium-neon  laser  is  used  in  the  laboratory  so 
the  single  mode  fiber  is  purchased  to  support  the  633  nm  wavelength.  Approximately  0.5 
m  of  single  mode  fiber  is  cleaved  for  the  sensor.  Tlie  competed  sensor  is  later  coupled 
to  the  experimental  setup  using  Norland  Lab  Splices.  The  fiber  that  provides  the  opposite 
mirrored  surface  is  a  multimode  fiber  because  of  its  lower  cost.  This  fiber  only  needs  to 
be  long  enough  to  insert  into  the  glass  tube  but  is  generally  cut  to  approximately  50  mm 
for  ease  of  handling. 

The  mirroring  is  done  with  sputtering  equipment  and  aluminum  is  used  to  create 
the  mirrors.  The  fibers  are  first  carefully  cleaved  and  cleaned  to  provide  a  flat  surface 
upon  which  to  sputter  the  aluminum.  Sputtering  times  may  vary  with  equipment  so  a  glass 
slide  is  first  placed  in  the  sputterer  to  determine  an  appropriate  time.  Reflectivity  should 
be  approximately  25%  on  the  single  mode  fiber  and  from  25%  to  100%  on  the  multimode 
fiber. 

After  the  fibers  were  coated,  the  glass  tube  was  prepared.  A  short  length  of  glass 
tubing,  approximately  25  mm,  was  cut  and  one  end  cleaved  fiat.  The  tubing  measured  135 
microns  on  the  inside  diameter  and  190  microns  on  the  outside  diameter.  Both  the 
prepared  tubing  and  one  of  the  coated  fibers  were  placed  into  a  jig  that  was  developed 
to  allow  alignment  of  the  fiber  and  the  tube  under  an  optical  microscope.  The  Jig  is 
shown  in  Fig.  2.4  and  consists  of  an  X-Y  positit>ner  with  Z  translation  and  a  mirror 
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holder/positioner  that  is  used  to  allow  the  fiber  and  the  tube  axes  to  be  made  parallel. 
This  jig  was  angled,  with  respect  to  the  horizontal,  to  allow  viewing  of  the  tube  end  while 
the  fiber  was  being  positioned  for  insertion.  This  open  tube  end  always  faces  up  and  is 
located  over  the  small  white  reflective  surface  to  provide  back-lighting. 

The  first  fiber  was  inserted  into  the  prepared  end  of  the  tube  and  bonded  using  a 
five  minute  epoxy.  The  fiber  was  inserted  so  as  to  achieve  the  desired  gage  length  and 
to  allow  approximately  10-40  microns  between  coated  fiber  ends  when  the  sensor  was 
completed.  After  the  epoxy  was  cured,  the  tube  was  cleaved  again  to  the  desired  gage 
length,  the  second  fiber  was  inserted  and  positioned  for  the  correct  gap,  and  the  epoxy 
was  applied  to  the  second  fiber.  Once  the  epoxy  has  cured,  an  accurate  measurement  of 
the  finished  gage  length  is  made  under  the  microscope,  and  the  sensor  is  ready  for  use. 


23  Calibration 

The  principle  of  operation  of  a  Mach-Zehnder  interferometric  fiber  optic  sensor 
was  mentioned  above.  The  axial  strain,  e,,  is  directly  proportional  to  the  number  of 
fringes,  N,  moving  past  a  fixed  point  and  can  be  given  as  (Butter  and  Hocker,  1978) 

t^=DN  (2.1) 

wfKic  D  is  the  proportionality  constant  and  can  be  written  as 


(2.2) 


u  :  t  -•  /,  is  the  vacuum  wavelength  of  the  optical  beam,  1  is  the  gage  length  cn'  the  sensor, 
n  is  the  refractive  index  of  the  fiber,  and  c  is  a  constant  which  depends  on  .strain  optic 
cau^ft'icients,  the  refractive  index,  and  Poisson’s  ratio  of  the  fiber.  It  is  possible  to 
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calculate  D  for  pure  silica,  however,  dopants  are  used  to  alter  the  refractive  index  of  the 
core  of  the  fiber  to  keep  the  light  propagating  in  the  fiber.  Hence,  the  fiber  has  to  be 
calibrated  to  determine  the  constant  D. 

The  Mach-Zehnder  sensor  was  calibrated  using  a  cantilever  beam  setup  as  shown 
in  Fig.  2.5.  The  fiber  sensor  was  mounted  to  a  cantilever  beam,  opposite  a  conventional 
strain  gage.  The  beam  was  displaced  and  a  plot  of  the  axial  strain  recorded  from  the  strain 
gage  vs.  the  fringes  was  made  for  a  particular  gage  length  of  the  fiber.  A  straight  line  was 
fit  through  the  data  points.  The  slope  of  this  straight  line  was  the  constant  D  as  shown 
in  Fig.  2.6.  The  strain  required  to  cause  one  fringe  to  move  across  a  fixed  point  is 


(2.3) 


inversely  proportional  to  the  gage  length  of  the  sensor.  In  other  words,  equation  (2.1)  can 
be  written  as  shown  in  equation  (2.3)  where  A  is  the  fringe  strain  sensitivity  constant  for 
the  sensor.  The  fringe  strain  sensitivity  can  be  defined  as  the  strain  required  to  cause  one 
fringe  to  move  across  a  given  point  for  a  unit  gage  length  of  the  sensor.  The 
experimentally  obtained  value  for  A  was  found  to  be  625  /zmm/mm  fringe  ‘  mm. 

Strain  is  obtained  from  the  Fabry-Perot  sensor  in  much  the  same  way  as  the  Mach- 
Zehnder  sensor.  However,  the  relationship  between  the  number  of  fringes  passing  and  the 
strain  is  only  dependant  on  the  wavelength  of  the  light  propagating  and  the  gage  length 
of  the  sensor.  Calibration  of  the  semsor  is  not  necessary.  The  relationship  between  the 
fringe  shift  and  the  axial  strain  is  explicitly  given  by  equation  (2.4).  Again,  N  is  the 
number  of  fringes,  X  is  the  wavelength  of  light,  and  L  is  the  measured  gage  length  of  the 
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sensor. 


(2.4) 


2.4  Evaluation  of  the  Fiber  Optic  Sensor 
2.4.1  Tension  Testing 

The  Mach-Zehnder  sensor  was  used  for  the  initial  evaluation  of  the  fiber  optic 
sensor.  A  single  mode  fiber  with  an  outer  diameter  of  80  microns  was  stripped  of  its 
plastic  jacket  in  the  region  that  was  to  be  bonded  to  the  specimen.  The  plastic  jacket  must 
be  stripped  to  ensure  proper  strain  transfer  between  the  specimen  and  the  optical  fiber. 
Also,  the  surface  of  the  specimen  must  be  treated  in  accordance  with  conventional  strain 


gage  techniques  prior  to  attachment.  Two  pieces  of  tape  were  placed  on  the  specimen 
leaving  the  desired  gage  length  exposed  between  the  tape  pieces.  These  tape  pieces  serve 
as  masks  so  that  only  the  desired  length  of  fiber  will  be  bonded.  The  tape  mask  was  then 
sliced  along  the  direction  that  the  fiber  will  run.  This  was  done  to  allow  the  tape  to  be 
removed  after  tlie  fiber  is  bonded.  Fig.  2.7  shows  the  specimen  with  the  tape  mask  in 
place.  The  fiber  was  then  laid  into  the  slice  and  bonded.  After  the  adhesive  cures,  the  tape 
is  removed. 

A  simple  tension  test  was  used  to  characterize  the  fiber  sensor  under  axial  loading 
conditions.  A  schematic  of  the  setup  for  the  Mach-Zehnder  interferometer  was  shown  in 
Fig.  2.1.  The  specimen  was  a  dog  bone  specimen  with  a  resistance  type  si  rain  gage 
bonded  to  one  side  of  the  specimen  and  the  fiber-optic  sensor  bonded  to  the  other  side, 
directly  opposite  the  strain  gage.  The  experiments  were  designed  for  axial  strain 
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measurements  only.  To  standardize  the  bonding  procedure,  both  sensors  were  bonded  with 
M-Bond  200,  a  cyanoacrylate  based  adhesive,  using  the  technique  described  in 
Measurements  Group  Bulletin  309 A.  Ambient  conditions  in  the  laboratory  are  maintained 
at  21  “C,  ±  3,  and  65%  relative  humidity,  ±  5. 

The  specimens  were  then  loaded  at  a  constant  strain  rate  and  the  strain  gage  signal 
was  recorded,  along  with  the  light  intensity  signal  from  the  photodiode,  with  a  Lecroy 
data  acquisition  system.  In  each  experiment,  loading  was  continued  until  failure  of  either 
the  strain  gage  or  the  fiber-optic  sensor. 

Another  series  of  experiments  was  conducted  in  which  the  specimen  was  loaded 
in  incremental  values  from  an  unloaded  condition  to  a  failure  condition.  The  specimen 
was  unloaded  and  removed  from  the  fixture  at  500  microstrain  increments  and  taken  to 
an  optical  microscope  for  inspection.  The  specimen  was  examined  and  photographs  were 
taken  of  any  significant  changes  in  the  adhesive  bond.  After  reconnecting  the  two  sensors 
to  the  data  acquisition  equipment,  the  specimen  was  strained  to  the  next  incremental  value 
and  the  procedure  was  repeated. 

Fig.  2.8  shows  typical  output  from  the  strain  gage  and  the  Mach-Zehnder  sensor 
for  the  axial  loading  experiments.  The  strain  data  from  both  the  strain  gage  and  the  fiber 
sensor  is  shown  in  Fig.  2.9. 

It  can  be  seen  that  the  Mach-Zehnder  sensor  was  linear  and  the  sensor  can  be  used 
for  strain  measurements.  However,  the  fiber-optic  sen.sor  was  found  to  fail  at  strain  levels 
far  below  the  levels  at  which  conventional  strain  gages  fail.  Failure  modes  could  be 
characterized  in  two  ways: 

•  the  fiber  fractures  and  the  sensing  signal  is  lost. 
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*  the  fiber  delaminates  with  matrix  cracking,  resulting  in  a  change  in 
effective  sensor  gage  length. 

Note  that  the  delamination  and  matrix  failure  is  related  to  the  bonding  of  the  fiber  and 
does  not  represent  a  shortcoming  of  the  glue.  A  separate  experiment  was  conducted  with 
only  a  thin  layer  of  glue  and  the  glue  did  not  show  cracking  or  failure  during  the 
experiment.  The  maximum  value  reached  during  this  experiment  was  2  percent  strain. 

The  data  presented  in  Fig.  2.9  shows  failure  by  fracture  of  the  fiber  which  resulted 
in  loss  of  the  sensing  signal  and  therefore,  loss  of  the  interference  fringes.  Both  modes 
of  failure  produced  unusable  data  at  values  between  1.2- 1.8  percent  strain.  Photographs 
taken  during  the  incremental  loading  are  shown  in  Fig.  2.10.  Fig.  2.10(a)  shows  the  fiber 
in  an  unstrained  state  and  as  it  appears  up  until  failure.  The  bright  line  running  down  the 
center  of  the  fiber  is  the  contact  between  the  fiber  and  the  specimen.  Note  that  the  glue 
does  not  cover  the  entire  fiber  but  rather  forms  a  bed  under  the  fiber. 

In  Fig.  2.10(b),  the  delaminated  region  shows  up  as  a  wide,  bright  band.  This  is 
caused  by  the  defracting  of  the  microscope  light  from  the  air  gap  that  opened  up  between 
the  fiber  and  the  specimen.  Under  actual  observation  the  band  appears  as  colored  fringes. 
Note  also  that  the  glue  matrix  shows  cracking.  These  cracks  will  always  curve  toward  the 
portion  of  the  fiber  that  is  unattached  and  unstrained,  and  therefore  show  that 
delamination  has  occurred  from  left  to  right  in  the  photograph.  The  right  side  of  the 
photograph  shows  the  portion  of  fiber  that  has  not  yet  delaminated. 

Fig.  2.10(c)  represents  an  experiment  in  which  the  fiber  fractured.  There  is  a 
central  crack  in  the  matrix  that  is  aligned  with  the  fracture  in  the  fiber.  All  other  cracks 
curve  toward  the  fracture  area.  This  matrix  cracking  shows  that  subsequent  delamination 


of  the  fiber  occurred  after  the  fracture  and  ran  outward  from  the  fiber  fracture  in  both 
directions. 

Fig.  2.11  shows  strain  data  obtained  from  the  incremented  experiment.  The  strain 
gage  data  was  plotted  using  a  linear  regression  curve  fit.  The  data  shows  that  the 
specimen  remained  linearly  elastic  throughout  the  experiment.  The  first  run,  represented 
by  the  circles,  and  the  second  run,  represented  by  the  inverted  triangles,  was  made  up  to 
1  and  1.5  percent  strain,  respectively.  Following  the  second  run,  observation  of  the  fiber 
sensor  showed  delamination,  see  Fig.  2.10. 

Analysis  of  the  raw  data  showed  the  effect  of  the  delamination  on  the  fiber  sensor 
output.  Fig.  2.12  shows  an  anomaly  that  occurred  at  approximately  13,200  microstrain, 
during  the  2nd  run.  The  fringes  from  the  1st  run,  where  no  anomaly  was  seen,  and  those 
that  precede  and  follow  the  anomaly  in  the  2nd  run,  are  regular  and  do  not  appear 
distorted. 

Since  the  strain  calculation,  given  in  equation  (2.1),  is  based  on  the  gage  length, 
the  observed  strain  will  be  affected  if  a  change  in  gage  length  is  not  accounted  for.  This 
affect  can  be  seen  in  Fig.  2.11,  runs  3  and  4.  Runs  3  and  4  were  made  up  to  1.6  and  1.7 
percent  strain,  respectively.  Fringes  from  both  of  these  runs  shov,  d  numerous  anomalies. 

From  these  observations,  it  is  believed  that  the  initial  delamination  occurred  at  the 
end  of  the  second  run  and  the  damage  grew  during  the  3rd  and  4th  runs.  This  damage 
growth  is  illustrated  by  the  difference  in  strain  observed  in  those  runs. 

From  the  initial  failure  experiments,  it  was  determined  that  the  fiber  sensor 
typically  failed  at  strains  around  1.5%.  The  first  strain  level  at  which  to  stop  was  set  at 
1%  strain,  or  10,000  microstrain.  In  all  experiments,  the  sensor  showed  no  significant 


changes  before  reaching  this  strain  level.  The  SOO  microstrain  increment  was  then 
followed  until  failure. 

From  these  characterization  experiments  it  was  determined  that  fiber  optic  sensors 
could  be  confidently  used  up  to  a  strain  level  of  1.S  percent. 

2.5  Experimental  Procedures  for  Contact  Strain  Measurements 

The  first  set  of  experiments  was  aimed  at  applying  fiber  optic  sensors  to  measure 
strains  due  to  a  static  load.  Disks  32  mm  in  diameter  were  machined  from  6  mm  thick 
plexiglass  sheets.  A  fiber-optic  Fabry- Perot  strain  sensor  was  bonded  to  the  surface  of  the 
disk,  close  to  the  point  of  the  applied  load.  The  disk  was  then  loaded  in  compression  and 
data  was  collected  from  a  PCB  load  cell  and  a  photodiode  using  a  Lecroy  data  acquisition 
system.  The  disk  was  loaded  to  approximately  1000  N. 

The  second  set  of  experiments  were  designed  to  measure  out  of  plane  strains  in 
a  disk  using  the  same  loading  configuration  as  that  used  with  the  Fabry-Perot,  however, 
a  Mach-Zehnder  sensor  was  used.  A  2  mm  hole  was  drilled  through  the  thickness  of  the 
disk.  The  jacket  on  the  fiber  was  removed  and  the  fiber  was  embedded  through  the 
thickness  using  a  chemical  called  Envirotex,  whose  mechanical  properties  are  very  similar 
to  those  of  plexiglass. 

Again,  a  piezoelectric  load  cell  placed  directly  below  the  disk  was  used  to  acquire 
the  loading  data.  From  theory,  strains  can  be  calculated  at  any  point  on  the  disk  if  the 
load  is  known.  As  mentioned,  the  fiber  was  embedded  on  the  loading  axis  and  close  to 
the  contact  point  so  that  Hertz  equations  could  be  utilized  for  the  analysis.  A  Ltcroy  data 
acquisition  system  was  used  to  record  both  the  fiber  optic  and  load  cell  signals. 


For  the  three  dimensional  experiments,  2S.4  mm  plexiglass  spheres  were  embedded 
with  a  fiber,  utilizing  the  Mach-Zehnder  sensor,  through  a  2  mm  hole.  The  experimental 
procedure  was  the  same  as  that  for  the  disks  with  the  embedded  fiber. 

2.6  Theory 

2.6.1  Two  dimensional  strains 

The  stress  field  equations  for  a  disk  are  given  in  (2.5),  (2.6),  and  (2.7)  (Smith  and 
Liu,  1953). 


(2.5) 
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The  terms  <|>[  and  <|»2  are  given  by  the  following 
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The  term  A  is  a  parameter  that  is  dependant  on  the  material  elastic  properties,  E  and  v, 
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and  the  radii  of  curvature  of  the  two  objects  in  contact,  and  is  given  in  equation  (2.8). 


*2  ) 


(2.8) 


Ute  friction  factor,  p,  is  taken  to  be  zero  because  there  is  no  tangential  load  and 
the  half  contact  width,  b,  is  found  using  the  load  data  and  given  by  equation  (2.9). 
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The  strain  along  the  loading  axis,  e^,  is  given  by  equation  (2.10).The  assumption 
of  plane  stress  conditions,  cyyy=0,  has  been  made  due  to  the  geometry  of  the  disks  and 
equation  (2.10)  reduces  to  equation  (2.11). 
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Equation  (2.11)  is  the  theoretical  strain  on  the  surface  of  the  disk  given  the 
geometrical  conditions  and  the  load. 

From  theory,  the  transverse  strain,  is  given  by  equation  (2.12). 


^yy  =  i(^rr  '  ^  a„)) 


(2.12) 


Again,  a  plane  stress  assumption  is  made  given  the  disk  geometry  and  (2.12)  is 
reduced  to  (2.13). 
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(2.13) 
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To  summarize,  if  the  toad  is  known  and  the  geometry  of  the  two  bodies  in  contact 
is  given,  the  surface  strain,  equation  (2.11),  and  the  transverse  strain,  equation  (2.13),  can 
be  calculated. 

If  the  fiber  is  far  away  from  the  diametrically  opposite  side,  the  effect  of  the  load 
from  the  opposite  side  is  relatively  small  and  this  has  been  shown  in  Fig.  2.13.  This  plots 
shows  the  percent  error  in  the  calculation  of  transverse  strain,  as  a  function  of  the 
normalized  position  of  the  sensor. 

2.6.2  Stress  in  a  sphere 

Johnson  (1985)  presents  equations  for  the  stresses  along  the  z-axis  in  a  sphere. 
These  equations  are 
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The  terms  in  these  expressions  are  given  in  the  following  equations. 
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The  effect  of  superposition  from  the  diametrically  opposite  load  was  neglected  for 


43 


• 

1=  1 . 1 
u  X 

(2.17) 

E-  ~X - ST 

(2.18) 

• 

_  _  3P 

(2.19) 

( CTy+o, ) ) 

(2.20) 

the  reasons  stated  above  in  the  theoretical  derivation  for  a  disk  under  diametral 
compression. 


2.7  Results  and  Discussion 

From  the  results,  it  is  seen  that  the  Mach-Zehnder  fiber  optic  sensor  can  be 
embedded  in  bodies  and  used  to  record  strains.  Fig.  2.14  and  2.15  are  theoretical  plots  of 
expected  strain  values  and  expected  number  of  fringes  for  varying  loads  and  varying  the 
position  of  the  embedded  sensor  on  the  z-axis.  Fig.  2.16  and  2.17  show  a  comparison  of 
the  experimental  and  theoretically  obtained  strain  values  for  the  disk  experiments.  Fig. 
2.18  shows  the  typical  outputs  of  the  fiber  optic  sensor  and  the  piezoelectric  load  cell. 
Fig.  2.19  is  a  plot  of  the  experimental  and  theoretical  half  contact  width  vs.  the  contact 
load.  A  consistent  deviation  in  the  slope  of  the  theoretical  and  experimental  plot  was 
observed  in  all  the  experiments  using  the  Mach-Zehnder  embedded  in  disks.  A  number 
of  reasons  for  these  errors  were  explored  and  are  listed  below: 

•  improper  measurement  of  the  gage  length  of  the  fiber.  Due  to  the  nature 
of  embedding,  the  measured  value  of  embedded  length  of  fiber,  or  in  other 
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words  the  gage  length,  might  not  be  accurate  because  the  effective  load 
carrying  length  of  the  fiber  is  dependent  upon  the  embedding  technique. 

•  the  equations  used  to  compute  the  strain  were  based  on  infinite  elastic  half 
spaces  with  a  Hertz  contact  area.  An  attempt  to  include  the  effect  of  the 
force  being  applied  at  the  diametrically  opposite  side  was  made  based  on 
the  work  in  Johnson(1985).  A  concentrated  force  was  assumed  for  the 
opposite  side  and  the  superposition  principle  was  applied.  The  theoretical 
plot  of  the  expected  error  is  shown  in  Fig.  2.13.  The  superposition  of  the 
opposite  force  did  not  contribute  significantly  to  the  strain  computations. 

•  the  Mach-Zehnder  sensor  yields  an  average  value  of  the  strain  over  the 
embedded  length  but  the  theoretical  equations  are  valid  for  a  point  on  the 
body.  The  development  of  an  integral  expression  for  strain  along  the 
embedded  length  might  produce  more  accurate  results.  Another  alternative 
would  be  to  use  a  "point"  sensor,  such  as  a  Fabry-Perot  sensor,  because  of 
its  ability  to  measure  strain  within  a  small  gage  length. 

•  the  presence  of  air  bubbles  in  the  filling  compound  might  have  altered  the 
effective  gage  length. 

Fig.  2.20  shows  the  data  obtained  from  the  Fabry-Perot  sensor  attached  to  a  disk 
under  diametral  compression.  Again,  the  theoretical  line  was  obtained  using  load  data 
from  the  load  cell  and  applying  Hertz  contact  theory  to  calculate  the  strain  at  the  sensor 
location.  The  plot  shows  agreement  with  the  theciretical  prediction  within  14%. 

Predicting  stress  and  strain  profiles  in  three  dimensional  bodies  is  an  order  of 
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magnitude  more  complicated.  For  the  three  dimensional  experiments,  a  sphere,  which  is 
the  most  simple  representation  of  a  3d  body,  was  used  as  the  specimen  and  a  fiber  was 
embedded  using  the  same  process  as  for  the  disks.  To  avoid  numerical  procedures,  strains 
were  computed  for  a  point  along  the  z-axis. 

As  was  the  case  in  disks,  a  very  significant  deviation  in  slope  was  observed 
between  the  theoretical  and  experimental  values  of  strain.  Once  again,  the  experimentally 
obtained  value  of  strain  is  the  average  value  of  strain  along  the  embedded  length  where 
as  the  theoretically  obtained  values  are  for  a  point  on  the  z-axis.  In  a  sphere,  the  strain 
increases  in  magnitude  from  the  free  surface  to  a  interior  point  on  the  z-axis(Johnson, 
1985).  Perhaps  an  arithmetic  mean  of  the  theoretical  strain  would  afford  better  comparison 
with  the  experimentally  obtained  value  of  average  strain  along  the  embedded  length. 

Fig.  2.21  shows  the  characteristic  fiber  optic  and  load  cell  output  for  a  sphere 
under  compression.  The  theoretical  equations  are  non  linear  and  a  slight  amount  of  non 
linearity  was  observed  in  the  strain  vs.  load  plots  as  seen  in  Fig.  2.22. 

2.8  Conclusion 

The  Mach-Zehnder  fiber-optic  sensor  measured  axial  strains  with  acceptable 
accuracy  and  these  strain  measurements  have  been  shown  to  be  linear  up  until  the  failure 
point  of  the  sensor.  This  failure  point,  using  the  technique  described  in  this  paper,  has 
been  shown  to  be  much  lower,  typically  1.2-1.8  percent  strain,  than  that  taken  to  be 
accepted  maximums  for  attached  resistance  type  strain  gages,  3-6  percent  strain.  Failure 
has  been  shown  to  be  either  due  to  fiber  fracture  or  delamination  and  glue  failure.  For 
applications  in  our  laboratory,  this  upper  limit  of  measurement  is  above  the  strains 
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encountered.  Experiments  were  also  conducted  in  which  the  entire  fiber  was  covered  with 
glue.  Covering  the  fiber  entirely  with  glue  produced  fiber  fracture  only.  The  glue  showed 
cracking  at  the  fiber  fiacture  point  and  the  failure  strains  were  comparable  with  the  values 
presented  in  this  paper. 

The  Mach-Zehnder  sensor  can  be  embedded  through  bodies  to  measure  strains. 
The  strains  predicted  by  the  sensor  are  an  average  over  the  gage  length  and  can  only 
approximate  the  actual  strain  value  at  a  unique  point.  The  technique  used  to  embed  the 
fibers  is  critical  for  accurate  prediction  of  strain  values.  The  experiments  using  the  32  mm 
diameter  disks  proved  that  transverse  strains  due  to  Poisson’s  ratio  effect  can  be 
approximately  predicted  by  these  sensors.  Sub  surface  strains  in  spherical  bodies  under 
diametral  compression  can  also  be  predicted  by  Mach-2!ehnder  fiber  optic  sensors. 
However,  the  value  of  the  strain  predicted  is  an  average  strain  over  the  embedded  length 
of  the  sensor.  This  strain  prediction  is  very  large  when  compared  to  the  theoretical  strain 
at  a  point  inside  of  the  sphere.  For  comparison  to  the  theory,  a  more  suitable  sensor 
would  be  the  Fabry- Perot  interferometric  sensor. 

Fabry-Perot  sensors  have  been  shown  to  be  applicable  for  determining  contact 
parameters.  The  sensor  was  within  14  percent  agreement  with  the  theoretical  values. 

The  mechanisms  of  failure  described  for  the  Mach-Zehnder  sensor  are  speculated 
to  be  indicative  of  those  for  the  Fabry-Perot  sensor.  However,  preliminary  work  has 
shown  that  there  is  a  dependance  of  the  failure  strain  on  both  the  diameter  and  the 
material  properties  of  the  glass  tube  used  in  constructing  the  Fabry-Perot. 

This  chapter  represents  preliminary  results  of  contact  strain  measurements  using 
fiber  optic  sensors.  The  work  will  be  contiunued  in  the  coming  year. 
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Figure  2.1.  Schematic  of  the  experimental  setup  for  the  Mach-2^hnder  sensor. 
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Figure  2.2.  Schematic  of  the  Fabry-Perot  sensor. 
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Figure  23.  Schematic  of  the  experimental  setup  for  the  Fabry-Perot  sensor. 


Figure  2.5.  Calibration  plot  for  the  Mach-Zebnder  sensor. 
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Fijjiui  e  2.6.  Schematic  of  cantilever  beam  for  calibration  of  the  Mach-Zehnder  sensor. 
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Figure  2.11.  Strain  gage  and  fiber  sensor  data  from  an  incremented  experiment  showing  the 
effect  of  fiber  delamination. 


Figure  2.12.  Output  from  both  the  a.)  strain  gage,  and  b.)  the  fiber  sensor  which  has  begun  to 
delaminate  (taken  from  2nd  run). 
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Figure  2.13.  Plot  of  percent  error  in  transverse  strain  as  a  function  of  sensor  position  and  disk 
diameter  due  to  the  far  field  load. 
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Figure  2.14.  Plot  of  theoretical  transverse  strain  vs.  contact  load  for  two  25.4  mm  diameter 
Homalite  100  disks  in  contact. 
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Figure  2.15.  Plot  of  theoretical  contact  load  vs.  fiber  optic  sensor  fringe  number  for  two  25.4 
mm  diameter  Homalite  100  disks  in  contact. 
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Figure  2.16.  Plot  of  theoretical  strain  and  strain  from  an  embedded  sensor  vs.  load  for  the 
disk  under  diametral  compression. 
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Figure  2.17.  Plot  of  theoretical  strain  and  experimental  strain  from  a  fiber  optic  sensor  for  a 
disk  under  diametral  compression. 
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Figure  2.18.  Typical  output  from  the  load  cell  and  the  fiber  optic  sensor  for  the  disk  under 
diametral  compression. 


Figure  2.19.  Plot  of  half  contact  width,  b,  vs.  load  for  a  Plexiglass*  disk  under  diametral 
compression. 
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Figure  2.20.  Typical  data  from  the  surface  attached  Fabry-Perot  sensor  as  compared  to 
theoretical  values 
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Figure  2.22.  Plot  of  radial  strain  vs.  compressive  load  for  a  25.4  mm  Plexiglass*  sphere  under  diametral  compression. 


CHAPTERS 


APPLICATION  OF  WHITE  LIGHT  SPECKLE  PHOTOGRAPHY 
TO  MEASURE  RIGID  BODY  DISPLACEMENTS 
AND  DEFORMATIONS 


This  chapter  presents  the  technique  of  white  light  speckle  photography  with  the 
intention  of  using  it  to  measure  displacements  around  contacting  particles.  Results  from 
preliminary  experiments  are  also  presented  and  discussed. 


3.1  Introduction 

3.1.1  History  of  Speckle  like  Phenomena 

When  a  beam  of  light  is  passed  through  an  isotropic  transparent  medium,  light  is 
scattered  from  each  point  along  its  path.  This  scattering  process  occurs  because  of 
inherent  optical  roughness,  in  other  words,  small  particles  in  suspension  inside  the 
material.  This  scattering  of  light  is  very  essential  for  the  working  of  the  speckle 
technique.  The  history  of  speckle  like  phenomena  predated  the  laser  but  engineering 
application  of  speckle  techniques  started  only  after  the  development  of  the  laser.  The  term 
‘speckle’  originated  as  a  result  of  holography  when  dark  spots  of  destructive  interference 
were  observed  in  holograms.  All  the  speckle  techniques  differ  only  in  the  nature  of  the 
speckle  pattern  generation.  A  double  exposure  of  the  speckle  pattern,  one  before 
deformation  and  one  after,  is  a  common  denominator  to  all  spi  kle  techniques.  Leendertz 
(1970)  established  the  foundation  principles  of  speckle  interferometry.  Archbold  and 
Ennos  (1972)  and  later  Duffy  (1974)  proposed  a  different  technique  known  as  speckle 
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photography  which  involved  the  use  of  a  single  divergent  laser  beam  and  recording 
double  exposed  photographs  of  the  speckle  patterns.  Barker  and  Foumey  (1976)  showed 
the  successful  use  of  this  technique  for  interior  displacement  measurements.  Of  late, 
speckle  photography  using  ordinary  whitelight  has  been  developed  by  Chiang  (1979, 
1980,  1990)  and  Asundi  (1982,  1983)  which  has  opened  a  wide  range  of  possibilities  for 
the  application  of  this  technique. 


3.2.2  Formation  of  ^Speckle’ 

When  light  from  a  coherent  source  such  as  a  laser  falls  on  an  optically  isotropic 
medium,  scattering  of  the  laser  light  occurs  due  to  the  optical  roughness  of  the  medium. 
This  scattered  light  travels  different  paths  and  hence  they  interfere  with  each  other 
forming  bright  and  dark  spots  of  constructive  and  destructive  interference  giving  a  speckle 
like  appearance.  These  spots  are  called  ‘speckles’.  This  scattered  light  speckle  pattern  acts 
like  a  random  grid  and  if  the  deformation  is  sufficient  to  cause  a  shift  of  the  speckle 
pattern  by  at  least  one  speckle  diameter,  the  double  exposed  speckle  pattern  will  behave 
like  a  complex  diffraction  grating.  Speckle  size  can  be  calculated  using  the  Rayleigh 
criterion.  This  gives  an  average  diameter  in  the  image  plane  of  s=1.2XF  (Barker  and 
Foumey,  1976)  where  F  is  the  aperture  ratio  of  the  lens.  Corresponding  speckle  size  in 
the  object  is  given  by 


1.2XF 

s= - 

m 


(3.1) 


where  m  is  the  magnification. 
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3.2  Speckle  Techniques 
3.2.1  Speckle  Interferometry 

When  two  optically  interfering  beams  of  light  are  used  to  generate  a  speckle 
pattern,  the  technique  is  called  speckle  interferometry.  The  two  sources  of  optically 
interfering  beams  can  be  a  laser  in  combination  with  an  appropriately  placed  mirror  which 
causes  a  phase  difference  necessary  for  interference. 


3.2.2  Speckle  Photography 

Speckle  photography  uses  a  single  beam  of  light  to  generate  a  speckle  pattern. 
When  a  laser  is  used  as  the  source  of  light,  a  speckle  pattern  is  generated  due  interference 
between  the  scattered  light  inside  the  material.  This  technique  is  called  laser  speckle 
photography.  When  a  speckle  pattern  is  artificially  created  and  ordinary  light  is  used,  it 
is  called  whitelight  speckle  photography. 


3J  Theory  of  Fringe  Formation  in  Speckle  Photography 

Speckles  are  formed  due  to  random  interference  of  coherent  light  from  different 
scattering  sites  within  the  propagating  medium.  The  intensity  of  an  individual  speckle  is 
independent  of  the  intensity  of  the  illuminating  light  and  depends  only  on  the  relative 
phase  of  the  light  from  each  scattering  site  within  a  speckle  cell.  Hence,  it  is  unimportant 
if  the  phase  of  the  illuminating  light  or  the  light  scattered  from  each  speckle  cell  changes 
between  exposures  due  to  material  birefringence  or  deformation  of  the  body.  The  double 
exposed  negative  of  the  speckle  pattern  or  specklegram,  as  it  is  called,  behaves  like  a 
complex  diffraction  grating.  When  a  laser  light  is  passed  through  this  grating,  light  is 
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scattered  into  a  diffraction  halo  with  intensity  varying  across  the  field.  The  diffraction 
pattern  is  due  to  the  product  of  the  diffraction  halo  of  a  single  speckle  pattern  and  the 
intensity  distribution  generated  by  two  point  sources  spaced  a  distance  apart  equal  to  the 
image  shift.  Sensitivity  is  governed  by  the  Rayleigh’s  criterion  and  the  maximum 
displacement  that  can  be  measured  is  governed  by  speckle  pattern  correlation.  Excessive 
deformation  moves  new  scattering  sites  into  and  old  scattering  sites  out  of  the 
illuminating  light.  Decorrelation  can  also  occur  due  to  a  displacement  gradient  normal  to 
illuminating  sheet. 

The  basic  optical  configuration  incorporating  a  single  beam  laser  speckle  method 
can  be  described  as  follows.  An  object  surface  is  illuminated  by  a  laser  beam.  The 
location  of  the  light  source  is  denoted  as  S  and  P  is  a  point  on  the  surface  of  the 
illuminated  body.  The  recorded  light  intensity  of  the  first  exposure  is  (Kobayashi,  1987) 


/,  -  B\x,y) 

If  the  body  is  deformed  between  exposures,  the  image  on  the  film  plane  will  be 
shifted  to  new  coordinates  x’  and  y\  The  light  amplitude  in  the  deformed  configuration 
may  be  written  as 


/,  -Bix',y') 

The  surface  displacement  vector  FP’  may  be  expressed  as 
Then  the  film  plane  displacement  vector  can  be  written  as 
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(3.4) 


PP'  -  ue  +  V  e  *  we 
p  *  p  y  p  t 

PP'flm  =  Vy] 

where  m  is  the  magnification  of  the  lens. 

Coordinates  of  the  film  plane  for  the  second  exposure  are  x’=  x  +  Uf  and 
+  Now  the  expression  for  the  total  intensity  is  expressed  as 

fpMai  =  *B\x  -y  Uf,y  *  (3-^ 

Interference  fringes  are  obtained  by  taking  optically  the  Fourier  transform  of  the 
amplitude  transmission  function  of  the  processed  photographic  negative.  The  amplitude 
transmission  function  g(x,y)  is  assumed  linear  in  the  ranges  of  interest  and  may  be 
expressed  as 

g{x,y)  =  a  +  bl^,  (3.7) 

where  a  and  b  are  constants. 

It  can  be  proved  that  the  intensity  in  the  transform  plane  can  be  written  as 


/(fa]^,o^)  =  C[(l-+cos(coMp)]  (3’8) 

where  C  is  a  constant  and  w  is  the  spatial  frequency  variable  and  is  denoted  as 


and 


to  = 


(3.9) 


(3.10) 


where  and  are  the  coordinates  in  the  transform  plane  and  K  =  2n/X.  Fringes  will  be 
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defined  when  =  0  i.e  when  1  +  cos({ou^  =  0. 


m 


For  a  small  region  of  illumination,  the  displacement  is  assumed  constant  and 
fringes  are  formed  when 

tou^  =  (2n  -  l)n  n  =  1,2,3...  (3*11) 

In  general,  only  the  first  fringe  is  observed  i.e  n  =  i,  therefore  written  as  in  equation 

(3.12),  where  d  is  the  distance  between  the  transfom  plane  and  the  specklegram. 

Xu.  +  Yv.  =  _  (3.12) 

»  R  ■»  /T  2 

To  find  the  component  of  displacement  along  F,,  we  take  X,  =  0  and  we  get 
=  )uil2Y^  and  similarly  =  Xdl2X^  Now  the  magnitude  of  total  displacement  is  given 
by 


(3.13) 


Substituting  values  of  the  components,  we  get  U  as  equation  (3.14), 


(3.14) 


where  S  can  be  calculated  as  shown  in  equation  (3.15). 


.S=- 


2XY 

5  S 


X^+Y^ 

*  s 


(3.15) 


In  case  of  a  magnification  while  recording  the  specklegram,  U  should  be  written  as 


u- 


tnkd 


(3.16) 
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3.4  Experimental  Procedure  and  Data  Analysis 

In  all  speckle  techniques,  a  double  exposure  of  the  speckle  pattern  is  made,  one 
before  deformation  and  one  after,  on  a  special  kind  of  film  called  the  litho  ortho  type. 
These  litho  films  are  a  go-nogo  kind  of  film  which  means  that  beyond  a  certain  intensity 
threshold,  they  will  print  as  black  and  below  the  threshold,  they  do  not  print  at  all.  In  the 
case  of  laser  speckle,  a  laser  source  is  used  to  generate  the  speckle  pattern.  In  the  case 
of  whitelight  speckle,  light  brom  a  Tungsten  filament  can  be  used  to  project  on  the 
artificial  speckle  pattern.  Fig.  3.1  shows  the  schematic  for  recording  a  whitelight 
specklegram.  These  artificial  speckle  patterns  can  be  a  rough  surface,  white  paint  sprayed 
on  a  black  base  or  a  coating  of  retroreflective  paint.  For  interior  measurements,  glass 
beads  can  be  mixed  while  casting  the  material.  A  double  exposure  of  this  speckle  pattern 
is  recorded,  one  before  deformation  and  one  after  deformation.  This  double  exposed 
speckle  pattern  is  called  a  specklegram.  Data  analysis  can  be  done  by  applying  either  the 
pointwise  filtering  technique  or  the  whole  field  filtering  technique.  When  a  laser  beam 
is  passed  through  any  point  inside  this  specklegram.  Young’s  fringes  are  seen.  This  is  the 
pointwise  filtering  technique.  The  spacing  between  the  fringes  is  related  to  the 
displacement  of  the  point  of  interest  as  shown  before  and  is  given  by  equation  (3.14).  Fig. 
3.2  shows  this  analysis  setup  which  is  called  the  pointwise  filtering  technique.  When  an 
optical  filter  with  a  pin  hole  is  used  for  analysis  of  the  specklegram,  whole  field 
displacement  information  can  be  obtained.  In  this  case,  isothetic  fringes  or  fringes  of 
constant  displacement  are  seen  and  the  analysis  is  termed  whole  field  filtering.  The  former 
technique,  though  tedious,  is  more  precise  and  can  be  automated.  A  computer  in 
combination  with  a  digital  filter  and  a  CCD  sensor  can  be  used  to  filter  the  noisy  speckle 
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data  and  give  accurate  displacement  values  at  different  points  in  the  model  (Chitsaz  and 
Moslehy,  1989). 

3.5  Results  and  Discussion 

Experiments  were  done  to  measure  the  surface  displacements  at  various  points  in 
a  plexiglass  disk  of  63.2  mm  diameter  subjected  to  rigid  body  displacement.  The  disk  was 
painted  black  and  then  sprayed  white  to  give  a  speckle  pattern  with  an  average  speckle 
size  of  around  100  /on.  The  disk  was  then  mounted  on  a  horizontal  translator  which  in 
turn  was  mounted  to  an  optical  bench.  2X90  W  Tungsten  lamps  were  used  to  project  light 
on  the  specimen.  Double  exposure  of  the  speckle  pattern  was  made  on  a  Kodalith  ortho 
type  3  film  using  a  4”X5"  field  view  camera.  In  between  the  double  exposures,  a  rigid 
body  displacement  of  400  fm  was  given  to  the  specimen.  The  film  was  developed  in 
Kodalith  Super-RT  developer.  The  specklegram  was  analyzed  using  the  pointwise  filtering 
technique.  A  1  mW  He-Ne  laser  was  used  for  the  analysis.  Fig.  3.3  shows  typical 
Young’s  fringes  and  Fig.  3.4  shows  the  displacements  calculated  at  various  points  on  the 
body.  The  calculated  displacements  compare  well  with  the  actual  displacement  of  400  /an. 
Most  of  the  measured  displacements  lie  within  10%  accuracy.  Accuracy  could  be 
improved  considerably  by  increasing  L,  the  distance  between  the  laser/specklegram  and 
the  viewing  screen,  especially  in  the  case  of  manual  measurement  of  fringe  spacing,  S. 

Experiments  were  also  done  using  25.4  mm  diameter  plexiglass  disks  under 
diametral  compression  using  the  crosshead  of  an  Instron  machine.  The  speckle  generation 
was  the  same,  but  the  fringe  formation  was  extremely  hazy  and  unclear  making  data 
acquisition  virtually  impossible. 
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Experiments  were  also  done  with  embedded  speckles.  A  three  dimensional  body 
was  made  out  of  a  polymer  resin  with  speckles  casted  inside  using  retroreflective  paint. 
The  experimental  procedure  was  the  same  but  the  fringes  again  were  very  hazy  and 
difficult  to  analyze. 

Speckle  techniques  are  powerful  means  for  the  non  destructive  evaluation  and 
testing  of  critical  components.  They  have  been  used  quite  extensively  in  the  turbine 
industry.  These  techniques  are  especially  suitable  for  measuring  interior  displacements 
which  is  a  limitation  of  most  techniques.  For  all  its  merits,  the  technique  has  its 
limitations  too.  Decorrelation  of  the  speckle  patterns  is  a  inherent  problem  of  laser 
generated  speckle  patterns.  This  can  be  avoided  by  using  whitelight  speckle  photography. 
White  light  speckle  technique  makes  it  possible  to  measure  displacements  at  any 
predetermined  point  of  the  body  while  avoiding  speckle  pattern  where  not  required.  The 
noisy  speckle  data,  which  is  a  inherent  defect  of  all  speckle  technique  can  be  analyzed 
accurately  by  using  a  digital  filter.  Moreover,  changing  the  white  light  source  to  an  ultra 
violet  source  might  improve  the  fringe  pattern  as  has  been  observed  by  some  researchers. 

3.6  Conclusion 

The  conclusions  that  were  drawn  from  the  speckle  experiments  can  be  summarized  as 
follows: 

•  Laser  speckle  photography  can  be  used  to  measure  displacements  on  the 
surface  of  a  body  although  decorrelation  of  speckle  patterns  is  a  problem. 

•  White  light  speckle  photography  is  able  to  measure  rigid  body 

displacements  on  the  surf  body  within  10%  accuracy. 


Noisy  speckle  patterns  were  observed  while  analyzing  speckle  patterns  of 
specimens  subjected  to  deforming  forces  and  embedded  speckles.  This 
problem  might  be  solved  by  a  digital  filter  and  by  changing  the  source  of 
whitelight. 


Currently,  experiments  are  being  designed  using  more  powerful  whitelight  sources 
and  digital  filters  to  analyze  the  speckle  data.  The  aim  is  to  embed  these  speckles  and 
gather  displacement  information  at  different  interior  points  in  the  body. 
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Fig.  3.1.  Experimental  setup  for  measuring  rigid  body  displacement  using  white  light  speckle 
photography. 


Fig.  3.2.  Data  analysis  setup  used  to  measure  rigid  body  displacements  in  a  plexiglass  disk. 


CHAPTER  4 


PRELIMINARY  STUDY  OF  THE  EFFECT  OF  CEMENTATION 

ON 

THE  DYNAMIC  LOAD  TRANSFER  PHENOMENON 

4.1  Introduction 

This  chapter  represents  a  preliminary  attempt  at  studying  the  effect  of  cementation 
on  the  wave  propagation  phenomenon.  Dynamic  wave  propagation  has  been  studied  in 
recent  years  for  its  importance  in  such  diverse  applications  as  dynamic  compaction  of 
powders,  shock  forming  of  modem  composites,  and  in  various  shock  absorbing  and 
isolating  applications.  These  studies  have  primarily  looked  at  the  particle  to  particle 
interaction  of  dry  contacts.  While  this  may  apply  to  the  fore-mentioned  applications,  the 
study  of  seismic  phenomenon  needs  to  account  for  other  parameters.  One  of  these 
parameters  is  cementation  between  particles. 

Predominantly  the  studies  to  date  have  been  based  on  the  classical  approach  to 
solutions  for  normal  compression  of  disks  and  spheres(Johnson,  1985).  This  approach 
addresses  the  problem  by  modeling  the  contact  as  that  between  two  elastic  bodies  where 
the  contact  zone  dimension  varies  with  varying  force.  If  cementation  is  considered,  the 
most  significant  difference  is  that  the  contact  zone  dimension  is  fixed  by  the  amount  of 
cementation  and  does  not  vary  during  the  application  of  force. 

In  nature  there  are  many  instances  of  cementation  beiw  een  particles.  Probably  the 
most  common  example  of  cementation  occurring  naturally  is  sedimentary  rock.  Man-made 
examples  include  such  things  as  concrete  and  asphalt.  By  studying  the  effect  of 


cementation  on  the  stress  wave,  we  can  gain  a  greater  insight  into  seismological  events 
and  structural  design.  Experimental  data  and  observations  are  needed  as  input  parameters 
for  numerical  simulations  of  wave  propagation  in  cemented  particulate  materials. 
Ultimately,  safer  design  and  construction  will  result. 

4.2  Theoretical  Analysis 

The  primary  analysis  technique  utilized  by  this  laboratory  is  an  overdeterministic 
least  square  method  to  determine  contact  parameters  from  digitized  photoelastic  fringe 
pattern  data.  The  theoretical  equations  are  based  on  Hertz  contact  theory.  The  equations 
for  the  stresses  were  given  in  Chapter  2  but  are  expanded  here  to  include  the  friction 
factor,  p. 
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The  terms  <j>,  and  <j>2  are  given  by  the  following 
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The  term  A  is  a  parameter  that  is  dependant  on  the  material  elastic  properties,  E  and  v, 
and  the  radii  of  curvature  of  the  two  objects  in  contact,  and  is  given  in  equation  (4.4). 


A  = 


1 

(l-v?  l-v?] 
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These  expressions  for  the  stress  field  are  combined  with  the  stress  optic  law  which 
is  shown  as  equation  (4.5), 


a. 


(4.5) 


where  N  is  the  fringe  number,  f,,  is  the  stress  optic  coefficient  for  the  material,  and  h  is 
the  thickness  of  the  specimen. 

Twenty  data  points  are  used  to  determine  the  two  unknowns,  b  and  p.  Once  the 
two  unknowns  have  been  found,  the  half  contact  width  and  the  friction  factor,  p,  are  used 
to  back  calculate  the  normal  and  tangential  contact  loads. 


43  Experimental  Procedures  and  Results 

Due  to  the  nature  of  the  cemented  contact,  static  experiments  were  conducted  to 
determine  whether  or  not  the  Hertz  contact  theory  would  be  valid  for  cemented  particles. 
As  mentioned  previously,  the  analysis  technique  utilized  by  this  laboratory  is  an 
overdeterministic  least  square  technique  in  conjunction  with  the  Newton-Raphson  method 


to  determine  contact  parameters  such  as  normal  load,  half  contact  width,  b,  friction  factor, 
p,  and  tangential  load,  based  on  Hertz  contact  theory.  The  theory  relates  the  contact  load 
to  the  contact  width.  Unlike  two  particles  in  contact  without  cementation  where  the 
contact  width  can  vary  with  the  load,  cemented  particles  have  a  fixed  contact  width. 

Static  photoelastic  experiments  were  conducted  using  two  circular  particles,  made 
of  Homalite-100,  that  were  cemented  together  with  a  cyanoacrylate  based  adhesive.  The 
disks  were  then  loaded  in  an  Instron  testing  machine  and  photographs  were  taken  through 
a  circular  polariscope  and  a  monochromatic  lens,  at  various  load  levels.  Fig.  4.1  shows 
typical  isochromatic  fringes  for  two  disks  with  a  cemented  contact,  under  diametral 
compression.  The  top  and  bottom  contacts  are  Hertz  contacts  and  the  disk  to  disk  contact 
is  the  cemented  contact. 

Looking  at  the  fringes  very  close  to  the  cemented  contact,  the  fringes  show  the 
effect  of  the  cementation.  Further  away  from  the  contact  the  fringes  are  very  similar  to 
the  Hertz  fringes.  It  is  this  similarity  that  was  hoped  to  be  utilized  in  later  dynamic 
experiments. 

The  disturbance  to  the  fringes  at  the  contact  point  seems  to  agree  well  with  the 
numerical  calculations  made  by  Dvorkin  et.  al.  (1991).  The  numerical  work  pi  Jicted  that 
for  very  stiff  cements  (stiffness  of  the  cement  is  expressed  with  respect  to  the  particle 
stifftiess)  the  profile  of  the  force  distribution  would  result  in  the  maximum  normal  forces 
at  the  edges  of  the  cemented  zone.  Dvorkin  also  calculated  the  normal  force  distribution 
for  soft  cements.  For  this  case  the  force  distribution  was  just  the  opposite,  the  maximum 
normal  force  would  fall  in  the  center  of  the  cemented  zone.  Fig.  4.2  shows  the  force 
distributions,  on  the  vertical  axis,  as  predicted  by  Dvorkin.  The  parameter,  m,  is  the  ratio 


of  the  cement  stif&iess  to  the  grain,  or  particle,  stiffness.  The  horizontal  axis  represents 
distance  along  the  particle  surface  from  the  center  of  the  contact  to  the  edge  of  the 
cemented  region,  ie.  0  -*  b. 

Static  photoelastic  experiments  were  designed  to  verify  these  predictions,  TTie 
disks  were  machined  of  PSM-1,  a  photoelastic  plastic,  which  has  a  low  stress  optic 
coefficient,  f^.  The  stress  optic  coefficient  determines  the  order  of  the  fringes  produced 
by  a  given  load.  The  higher  the  stress  optic  coefficient,  the  lower  the  fringe  order  seen 
for  a  given  load.  (The  lower  the  stress  optic  coefficient,  the  higher  the  fringe  order  for 
the  same  load.)  The  PSM-1  was  chosen  so  that  under  very  small  loads,  an  appreciable 
•  'amber  of  fringes  would  still  be  seen. 

The  experiment  called  for  a  uniform  loading  conditions  on  both  sides  of  a  disk  so 
three  disks  were  used.  The  center  disk  would  experience  a  distributed  load  on  both  sides, 
while  the  two  outer  disks  experience  a  point  loading  on  one  side  and  a  distributed  loading 
on  the  other.  To  fabricate  the  specimens,  the  cement  materials  were  either  cast  from  a 
liquid  into  the  gap  between  the  disks,  in  the  case  of  the  soft  cement,  or  machined  from 
a  solid  piece  and  ’glued’  into  the  gap,  for  the  stiff  cement.  The  procedures  will  be 
described  below. 

The  soft  cements  were  Hexcel  urethanes  of  varying  stiffness  that  were  mixed  from 
two  liquid  components  and  poured  between  the  two  disks.  The  disks  were  adhered,  face 
down,  to  a  sheet  of  plastic  with  the  desired  gap  between  them  and  the  width  of  the 
cementation  zone  was  fixed  by  placing  two  blocks  against  the  disk  edges.  The  area  was 
sealed  with  silicone  caulk  and  the  urethane  was  cast  into  the  resulting  cavity. 

For  the  stiff  cements,  a  piece  of  aluminum  was  machined  into  the  same  shape  as 


the  cast  urethane.  The  aluminum  ’cement’  was  then  bonded  to  the  PSM-1  disks  using  a 
photoelastic  adhesive.  The  adhesive  cures  to  the  same  mechanical  properties  as  the  PSM-1 
and  is  also  optically  clear. 

Fig.  4.3  and  4.4  show  the  typical  isochromatic  fringes  obtained  for  static 
experiments  involving  large  cementation  zones  and  for  different  cement  stiffness.  Fig.  4.3 
shows  the  fringes  that  result  when  the  cement  stiffness  is  very  low  compared  to  the  disk 
stiffness.  Notice  that  the  fringes  develop  from  the  center  of  the  cemented  region.  Just  the 
opposite  is  seen  for  the  case  of  the  very  stiff  cement,  Fig.  4.4.  The  fringes  develop  from 
the  edges  of  the  cemented  region.  Both  cases  confirm  the  numerical  prediction  made  by 
Dvorkin.  These  fringes,  however,  are  very  different  from  Hertz  fringes.  It  was  seen  that 
for  cementation  zones  that  were  large,  20%  of  the  particle  diameter,  and  3  mm  wide,  the 
fringe  patterns  deviated  more  from  the  Hertz  patterns  and  Hertz  analysis  would  not  be 
appropriate. 

Ihe  photographs  from  the  static  experiments  using  the  Homalite-100  disks  were 
analyzed  for  both  the  known  Hertz  contact  and  the  cemented  contact.  These  experiments 
involved  a  very  thin  cement  layer  and  also  the  cemented  region  was  not  wide  compared 
to  the  particle  diameter  (the  cemented  region  was  10%  of  the  particle  diameter).  It  was 
hoped  that  the  Hertz  equations  might  be  valid  for  some  region  and  analysis  could  be 
performed  using  the  previously  ment:  .;  .d  scheme. 

Typical  plots  showing  the  results  of  the  analysis  are  shown  in  Fig.  4.5.  The  plots 
show  the  experimentally  determined  load,  normalized  by  the  Instron  load  ceil  data,  vs.  the 
normalized  distance  from  the  contact,  i  he  normalizing  parameter  for  the  distance  from 
the  contact  was  the  half  contact  width,  h.  as  taken  from  the  Hertz  contact.  The  data  from 
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the  Hertz  contact  is  shown  as  solid  circles  and  the  dashed  lines  are  the  allowable 
experimental  error.  The  Hertz  data  shows  acceptable  values  within  a  region  from  about 
2b  to  8b.  This  agrees  well  with  previous  work  (Shukla  and  Damania,  1988). 

The  data  from  the  cemented  contact  shows  some  agreement  with  the  Hertz  data 
but  the  band  at  which  acceptable  values  can  be  obtained  is  very  narrow.  At  the  highest 
load,  2220  N,  the  data  falls  within  the  acceptable  error  between  approximately  Sb  and  7b. 
For  lower  loads,  however,  the  acceptable  region  shifts  to  higher  multiples  of  the  half 
contact  width,  b. 

The  narrowness  of  the  band,  as  compared  to  the  Hertz  data,  also  raises  some 
concern.  If  the  valid  region  for  analysis  is  too  narrow  and  shifts  with  the  applied  load,  the 
usefulness  of  the  analysis  technique  may  be  in  question.  However,  it  should  be  noted  that 
the  data  presented  does  not  represent  a  complete  study  and  further  work  will  be  done  in 
this  area.  An  acceptable  region  may  be  found  that  will  allow  analysis  of  the  isochromatic 
fringes. 

Following  the  static  experiments,  dynamic  experiments  were  begun.  A  single  chain 
of  cemented  particles  was  placed  in  the  camera  and  an  explosive  charge  was  detonated 
to  initiate  the  stress  wave.  Fig.  4.6  shows  typical  photographs  of  the  disk  chain  as  the 
dynamic  wave  propagates.  From  the  photographs,  the  pulse  length  can  be  directly 
measured  and  the  group  wave  velocity  can  be  determined.  The  pulse  length  is  marked  in 
the  figure. 

The  pulse  length,  shown  in  Fig.  4.6,  can  be  seen  to  be  approximately  four  and  one 
half  particle  diameters.  This  differs  from  findings  repxirted  in  earlier  work  (for  example; 
Shukla  and  Damania,  1987)  where  the  pulse  length  consistently  stabilized  at  four  particle 


diameters. 


The  group  wave  velocity  was  determined  by  measuring  the  position  of  the  wave 
front  at  various  times  and  plotting  the  distance  propagated  vs.  time.  Fig.  4.7  shows  the 
plot  for  a  number  of  single  chain  experiments  of  cemented  particles  in  a  vertical 
arrangement.  The  average  velocity  was  determined  to  be  1240  m/s.  This  is  higher  than 
the  reported  velocity  for  a  single  chain  of  particles  without  cementation  by  approximately 
20%.  The  velocity  for  particles  without  cementation  was  found  to  be  approximately  1050 


Fig.  4.8  shows  the  normalized  contact  load  vs.  time  for  the  single  chain  of 
cemented  particles.  Hertz  analysis  was  used  to  acquire  the  load  data  presented  and 
identical  methodology  was  used  for  analyzing  each  of  the  contacts.  The  normalization 
removes  the  quantitative  aspects  and  therefore  the  plot  shows,  qualitatively,  the  effect  of 
the  cementation  on  the  contact  load.  The  plot  shows  that  across  four  particle  contacts, 
there  is  no  significant  attenuation  of  the  contact  load.  This  matches  the  results  for  a  single 


chain  of  particles  without  cementation. 

Experiments  were  also  conducted  with  the  particles  staggered  from  the  vertical, 
by  some  angle,  a.  Fig.  4.9  shows  the  particle  arrangement  for  the  staggered  configuration. 
To  date,  only  two  different  angles,  a  =  15°  and  a  =  40°,  have  been  used  in  the 
experiments.  Fig.  4.10  and  4.11  show  the  typical  isochromatic  fringes  obtained  using 


these  particle  arrangements. 

The  wave  speeds  have  been  determined  and  are  shown  in  Fig.  4.12.  These  velocity 
measurements  repre.sent  the  vertical  velocity  of  the  wave  front  as  it  propagates  in  the 
chain.  Since  the  data  only  represents  a  few  experiments,  the  reported  velocities  are  still 
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considered  estimates.  The  trend  however  would  lead  one  to  believe  that  the  greater  the 
angle,  the  lower  the  wave  velocity.  It  is  difficult  to  discuss  pulse  length  in  terms  of 
particle  diameters  but  it  is  easy  to  see  that  the  pulse  is  loading  approximately  6  particles 
in  Fig.  4.10,  the  15°  angle  case,  and  5  particles  in  Fig.  4.11,  the  40°  angle  case. 

Notice  that  there  was  energy  still  propagating  in  the  particles  (this  is  most 
noticeable  for  the  40°  angle  case)  after  the  main  pulse  had  passed.  This  energy  is  due  to 
reflections  taking  place  in  the  particles  and  presumably  would  lead  to  greater  attenuation 
of  the  pulse.  Evidence  of  this  can  be  seen  in  Fig.  4.11.  The  magnitude  of  the  fringe  orders 
seen  in  the  photographs  decreases  dramatically  as  the  loading  pulse  travels  down  the 
chain.  As  complete  confidence  in  the  current  analysis  technique  has  not  yet  been 
established,  no  load  data  will  be  presented  at  this  time. 

Another  set  of  experiments,  still  in  the  preliminary  stages,  is  the  investigation  of 
stress  wave  propagation  at  angles  of  90°  and  greater  from  the  incident  wave  direction. 
The  particles  in  these  experiments  were  cemented  together.  Fig.  4.13  shows  typical 
isochromatic  fringes  from  such  an  experiment.  Note  that  in  previous  work  by  members 
of  this  laboratory  (Shukla  et.  al.,  1988)  the  wave  would  not  propagate  at  angles  of  90°  or 
greater  using  uncemented  particles. 

Fig.  4.13  shows  that  energy  does  indeed  transfer  to  the  particle  chain  that  is  at  a 
right  angle  to  the  incident  pulse.  A  careful  look  at  these  photographs  reveals  that  damage 
growth  has  begun  in  the  first  disk  in  the  right  angle  chain.  Tlie  stress  on  the  contact  is 
almost  entirely  shear  stress  and  ilu*  damage  originates  at  the  edge  of  the  cementation 
zone.  The  origin  of  »he  damage  ci)incides  well  with  the  numerical  predictions  made  by 
Dvorkin  (1991)  with  regard  to  the  maximum  stress  location  for  a  cemented  contact. 
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Velcx:ity  measurements  were  also  made  for  the  pulse  propagating  in  the  near 
vertical  straight  chain  and  it  was  found  that  the  ratio  of  cementation  width  to  particle 
radius,  w/r,  also  has  an  effect  on  the  velocity  of  the  pulse.  The  velocity  was  found  to  be 
higher  in  this  straight  chain  (w/r  =  0.12),  approximately  1330  m/sec,  than  that  of  the 
previous  straight  chain  (w/r  =  0.10),  1240  m/sec.  The  difference  between  the  two  chains 
was  that  the  ratio,  w/r,  was  higher  in  the  chain  with  the  higher  velocity.  Further 
experimentation  will  be  done  to  explore  this  phenomenon. 

4.4  Conclusion 

Cementation  proved  to  have  an  impact  on  the  nature  of  wave  propagation 
phenomenon.  The  velocity  of  the  group  wave  as  compared  to  uncemented  contacts 
increased  by  approximately  20%,  for  vertical  chains  of  particles.  This  preliminary  work 
shows  that  although  the  group  wave  velocity  and  the  pulse  length  seem  to  be  affected,  the 
contact  load  does  not  seem  to  show  any  appreciable  attenuation  in  the  straight  particle 
chains. 

The  effect  seen  on  the  group  wave  velocity  would  seem  to  be  consistent  with  the 
physical  changes  made  to  the  medium.  The  cementation  is  bringing  the  particulate 
medium  closer  to  a  continuous  medium  and  the  wave  velocity  is  increasing  accordingly. 

The  change  in  the  pulse  length,  when  compared  to  uncemented  particles,  needs  to 
be  explored  further.  Long  chain  experiments  with  strain  gages  will  be  conducted  to 
determine  the  entire  effect. 

Hie  effect  of  contact  angle,  while  not  complete,  would  seem  to  indicate  that  the 
angle  affects  the  velocity,  the  pulse  length,  and  the  load  attenuation.  One  possible 
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explanation  for  the  changes  to  the  velocity  and  attenuation  is  that  the  pulse  undergoes 
more  scattering  and  dispersion  upon  encountering  particle  boundaries. 
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Figure  4.2.  Plots  of  normal  stress  along  the  cement  region  as  a  function  of  the  ratio  of  cement 
stiffness  to  particle  stiffness,  (from  Dvorkin,  1991) 
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Figure  4.4.  Typical  isochromatic  fringes  for 
’stiff’  cementation. 
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Figure  4.5.  Plots  of  normalized  load  vs.  normalized  distance  from  the  contact  for  static 
experiments. 
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Figure  4.6.  Typical  isochromatic  fringes  obtained  from  a  single,  verticle  chain  of  cemented 
particles. 
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Figure  4.7.  Distance  vs.  time  plot  for  dciermination  of  group  wave  velocity. 
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Figure  4.8.  Plot  of  normalized  contact  load  vs.  time  for  the  single,  verticle  chain  of  particles. 
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Figure  4.9.  Sdunialic  the  disk  chain  used  for  the  cemented,  angled  contact  experiments. 
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Figure  4.10.  Typical  isochromatic  fringes  obtained  from  a  single,  15“  angle  chain  of  particles 


CHAPTER  5 

HYBRID  THEORETICAL/EXPERIMENTAL 
DETERMINATION  OF  CONTACT  LOADINGS  BETWEEN 
CEMENTED  DISK  PARTICLES 


5.1  Introduction 

The  determination  of  the  correct  contact  Jaw  is  a  Jcey  problem  in 
understanding  the  dynamic  response  and  wave  propagational  behaviors  of  granular 
materials.  In  this  regard,  our  previous  research  (Shukla,  et.al.  1985,  1987,  1988) 
using  dynamic  photoelasticity  has  been  very  useful  to  determine  the  transient  contact 
stress  field  inside  model  disk  particles,  and  has  provided  fundamental  information  to 
help  develop  appropriate  interparticle  contact  laws  that  have  been  used  in 
computational  models  (Sadd,  et.al.  1993).  To  date  our  studies  have  focused  on 
particulate  systems  with  dry  and  uncemented  contacts.  Interest  has  now  developed 
into  the  dynamic  behavior  of  cemented  particulate  materials,  and  the  purpose  of  this 
chapter  is  to  outline  our  initial  research  in  establishing  a  hybrid 
experimental/theoretical  scheme  to  determine  the  dynamic  load  transfer  through  such 
cemented  systems. 

Our  photoelastic  observations  indicate  that  the  local  contact  stress 
distributions  between  cemented  disk  particles  are  different  from  that  of  Hertzian 
contacts.  This  fact  follows  from  photoelastic  experiments,  where  it  is  found  that  the 
distribution  of  local  fringe  patterns  of  cemented  and  dry  disks  are  different. 
Qualitatively,  for  high  modulus  cements,  the  maximum  order  fringe  grows  from  the 
two  ends  of  the  contacting  surface,  while  soft  cements  produce  an  evolving  maximum 
order  fringe  interiorly  located  near  the  center  of  the  contacting  surface,  see  Figs, 
5.1(a)  and  (b).  For  the  dry  case  with  no  cementation  (presumed  to  be  Hertzian 
contact),  the  maximum  order  fringe  grows  from  the  center  of  the  contacting  surface 
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as  shown  in  Fig.  5. 1(c).  Thus  the  local  stress  distributions  close  to  the  contact  surface 
are  clearly  different  between  these  cases.  For  Hertzian  contact,  an  elliptical  contact 
stress  distribution  is  assumed,  and  the  evaluation  of  the  stress  field  within  the 
cemented  disks  can  be  obtained  from  classical  elasticity  theory,  (e.g.  Timoshenko,  and 
Goodier,  1970).  Results  of  such  a  calculation  (see  Fig.  5.2)  indicate  that  the 
predicted  maximum  shear  stresses  using  Hertz  theory  are  much  higher  than  the 
photoelastic  data  for  the  case  of  disks  cemented  by  a  stiff  material  (super  glue).  As 
shown  in  the  figure,  the  difference  reaches  50%  at  the  maximum  shear  stress  point  for 
the  contact  loading  of  3 121b.(1389N). 

Dvorkin  et  al,  (1991)  proposed  an  analytical  model  to  predict  the  influence  of 
the  cementation  on  the  interparticle  load  transfer  between  .  Using  elasticity  theory 
along  with  particular  assumptions  of  the  displacement  distribution  in  the  cement  layer, 
an  integral  equation  for  the  stress  distribution  on  the  contacting  surface  was 
determined.  This  solution  was  not  in  the  explicit  form  and  required  considerable 
numerical  calculations  to  determine  the  resultant  forces  between  the  cemented  disks. 
Dvorkin  et  al.  determined  the  stress  distributions  on  the  contacting  surface  for  a 
special  case  in  which  the  two  cemented  disks  touch  only  at  the  center  point  of  the 
contact  line.  These  results  show  that  the  contact  stress  distribution  profile  changed 
significantly  when  the  cement  material  modulus  varied  from  soft  to  hard.  When  the 
cement  layer  is  softer  than  the  disk  material,  the  maximum  stress  in  the  contact 
distribution  lies  in  the  center  of  the  contact  zone.  On  the  other  hand,  if  the  cement 
material  has  a  higher  modulus  than  the  disk  particle,  the  maximum  stress  occurs  at  the 
edges  of  the  contact  area.  Unfortunately,  this  scheme  would  not  be  very  useful  to  be 
used  with  photoelastic  data  to  develop  a  hybrid  method  to  process  dynamic  fringe 
data  to  calculate  interparticle  load  transfer. 

In  an  effort  to  develop  such  a  hybrid  scheme,  and  guided  by  Dvorkin's  results, 
various  trial  normal  stress  distributions  on  the  contacting  surface  have  been  used  to 
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calculate  the  maximum  shearing  stresses  within  contacting  disks.  The  assumed 
contact  loading  distributions  included  convex  parabolic,  elliptical,  uniform  and 
concave  parabolic  profiles.  Comparing  the  results  of  these  trial  contact  loadings  with 
photoelastic  experimental  data,  it  was  found  that  the  concave  parabolic  profile 
produced  an  acceptable  comparison  with  the  particular  data  set,  see  Fig.  S.2.  This 
contact  stress  distribution  was  then  used  to  compute  the  interior  stress  field  within  the 
disk,  and  this  information  was  used  to  generate  contour  plots  of  the  maximum  shear 
stress  (identical  to  photoelastic  fringe  patterns)  as  shown  in  Fig. 5. 3.  It  can  be  seen 
that  the  numerical  results  closely  match  the  experimental  fringe  patterns  shown  in  Fig. 
5.2(b).  These  preliminary  trial  results  indicate  that  particular  contact  stress 
distributions  on  the  cemented  contacting  surface  could  be  used  to  determine  the 
correct  stress  field  within  the  disk  particle.  Thus  an  inverse  problem  solution  strategy 
could  be  established  whereby  the  resultant  force  between  pairs  of  cemented  disks 
could  be  determined  from  photoelastic  experimental  data.  The  present  study  uses 
such  an  approach  to  calculate  the  contact  force  between  two  cemented  circular  disks 
using  a  polynomial  family  of  stress  distributions  on  the  contacting  surface.  The 
advantage  of  this  method  is  computational  efficiency  and  the  fact  that  characteristics 
of  the  cementation  need  not  be  known.  The  information  needed  are  the  positions  and 
the  fringe  orders  of  the  photoelastic  fringes  within  a  contacting  disk  particle,  the 
material  and  geometrical  properties  of  the  disk,  and  the  cemented  contact  width. 
Several  examples  are  given  to  verify  the  accuracy  of  this  hybrid 
theoretical/experimental  approach. 

5.2  Basic  Theory  and  Methodology 
5. 2. 1  Classic  Theory 

From  classical  elasticity  theor>  (Timoshenko  and  Goodier,  1970),  the  stress 
field  for  the  plane  problem  of  a  circular  disk  loaded  with  a  pair  of  equal  and  opposite 
forces  as  shown  in  Fig.  5,4,  is  given  by 


96 


a.=  j.sin^e^ose.^sini^^s92|^^cP3y 


,  _  2/?/sin0icos^0i_sin02COs202 
^  -  -Tr\ - n - n - 


(5.1) 


where  the  coordinates  are  indicated  in  Fig.  S.4. 

The  solution  for  a  distributed  boundary  loading  y?(p^,  as  shown  in  figure  5.5, 
can  be  obtained  by  integrating  the  fields  specified  by  equations  (5. 1)  to  get 
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o.=  j{-24^{5^+S^}  +  =^cos<p}«d.p 

Xxy=  j  ^/(^)  |sin0icos^0i_sin02cos^02|;^^p 

-<X 

The  maximum  shear  stress  is  given  by  the  well-known  result 


(5.2) 


T 


(5.3) 


Equations  (5.2)  and  (5.3)  indicate  that  can  be  calculated  if  the  boundary 
loading  distribution  function  /((p)  and  angle  a  are  known.  These  parameters  are 
directly  related  to  the  stresses  transmitted  to  the  disk  particle  from  the  cement  layer. 
5.2.2  Contact  Stress  Distribution 

Dvorkin  et  al.  (1991)  studied  the  stress  distributions  at  the  surface  of  two 
cemented  circular  grains  for  the  case  with  different  ratios  of  cement  to  grain  moduli. 
When  the  cement  is  soft  relative  to  the  particle,  the  maximum  shear  stress  takes  place 
in  the  centei  ot'  the  contacting  surface,  but  when  the  cement  is  stiff  the  maximum 
shear  stresses  occur  at  the  ends  of  the  contact  surface.  Based  on  these  results,  it  is 
assumed  that  all  contact  stress  distributions  could  be  represented  by  a  family  of 
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polynomials  containing  several  parameters  to  provide  the  appropriate  shape  for  the 
cement/particle  moduli  under  study.  Such  a  stress  distribution  polynomial  can  be 
written  in  the  following  form 


/((p)  =  P[l  +  (m-l) 


M 


] 


(5.4) 


where  P  is  the  stress  at  the  center  of  the  contacting  surface,  m  is  the  ratio  of  the 
center  to  the  edge  stresses,  and  n  is  the  power  of  the  polynomials. 

If  the  parameters  n  and  a  are  fixed,  and  the  resultant  force  between  the  disks 
is  held  constant,  the  parameter  m  will  dominate  the  appropriate  shape  of  the  curve 
(see  Fig.  5.6).  If  m  =  1,  the  distribution  is  uniform;  for  /n  >  1,  the  curve  is  concave; 
and  m<  1,  the  curve  is  convex.  Therefore,  the  loading  shapes  available  with  this /(((>)- 
family  approximately  coincide  with  the  results  obtained  by  Dvorkin  for  various 
cement/particle  moduli.  The  stresses  within  the  disk  particle  that  are  determined  by 
using  the  above  contact  stress  distributions  are  very  close  to  the  experimental  data.  If 
all  the  parameters  in  the  polynomial  are  given,  the  resultant  force  on  the  surface  is 
determined  by  integrating  the  polynomial  over  the  surface,  i.e. 


F,  =  2j/(<|>)«</(p  = 
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2PaR 


m  +  n 
n  +  \ 


5. 2. 3  Numerical  Scheme 


(5.5) 


Substituting  the  equation  (5.4)  into  (5.2)  and  (5.3),  produces  the  general  form 

=  9I(P,/n,//)  (5.6) 

which  is  an  implicit  expression  in  terms  of  p,  m  and  n. 

According  to  photoelasticity  theory,  an  isochromatic  fringe  pattern 
corresponds  to  a  line  of  constant  maximum  shear  stress  which  is  specified  by 
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where  /„  is  the  material  fringe  value,  N  is  the  fringe  order  number,  and  h  is  the 
thickness  of  the  disk  particle. 

Using  these  results,  an  error  function  g  can  be  defined  as 

f  -N 

Sk  =  - -  (5.8) 

where  the  subscript  k  refers  to  ^-th  data  point.  It  is  desired  for  this  error  function  to 
be  reduced  to  zero;  however,  applying  this  equation  to  actual  data  and  using 
approximate  analytical  models  for  the  calculated  stress  distributions  will  yield  non¬ 
zero  values  for  g^.  In  order  to  reduce  this  error  function  to  acceptable  levels, 
numerical  schemes  have  been  developed.  It  has  been  established  (Shukla  and  Nigam, 
1985;  McCalla,  1967;  Singer,  1964)  that  by  using  Gauss-Legendre  quadrature,  along 
with  a  Newton-Raphson  iteration  method,  a  numerical  procedure  can  be  constructed 
to  reduce  g^  . 

Applying  these  techniques  to  the  current  problem  involving  equations  (5.6)  - 
(5.8),  it  is  first  noted  that  initial  estimates  of  the  parameters  P,  m  and  n  usually  will 
not  make  gi=0,  and  an  iterative  scheme  is  required.  Equation  (5.8)  can  be  written  in 
the  form  of  a  first  order  Taylor  series  expansion  of  g^  ,  i.e. 


where  the  subscript  /  refers  to  an  /-th  iteration. 

Applying  this  result  to  numerous  photoelastic  data  points  produces  an  over¬ 
determined  system  of  equations.  Solving  this  over-determined  set  produces 
corrections  A/^,  Aw;,  and  A//^  which  are  used  to  improve  the  parameters  through  the 
relations 

P,.>  =  P.  +  ^  ^  =  ",  +  A",  (5  i0) 

for  the  next  iteration  step.  This  procedure  continues  until  the  error  difference  between 
iterations  is  smaller  than  a  specified  controlling  value. 
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5.3  Numerical  Results 

A  computer  code  implementing  the  above  procedures  has  been  developed,  and 
a  preliminary  numerical  investigation  has  been  conducted.  Results  from  the  code 
indicated  that  the  determination  of  the  resultant  cementation  forces  was  not  very 
sensitive  to  the  n  parameter,  and  in  fact  values  of  n  from  2  to  6  do  not  produce  any 
significant  changes  in  the  results.  Because  of  this  fact,  the  n  parameter  was  set  at  a 
constant  value  of  4,  and  only  the  parameter  p  and  m  were  varied  in  the  iterative 
scheme. 

Data  from  several  sets  of  photoelastic  experiments  were  used  to  verify  the 
hybrid  approach  previously  described.  The  cement  materials  that  were  used  in  the 
experiments  included  urethane,  super-glue,  plexiglas  and  aluminum,  which  provided  a 
variety  of  stiffnesses.  All  the  model  disks  were  made  of  polycarbonate.  The  data 
collection  was  completed  by  a  digitizer  that  gives  the  position  and  fnnge  order  of 
each  photographic  data  point.  The  choice  of  the  data  points  will  affect  the  accuracy  of 
the  results.  Generally  one  should  take  more  data  near  the  maximum  shear  stress 
points  because  the  variation  of  these  fringe  positions  with  load  is  more  sensitive  than 
those  of  other  fringes.  For  the  urethane  case  (soft  cement),  the  fringes  were  very 
thick,  see  Fig.  5.1(a),  and  it  was  difficult  to  determine  the  exact  fnnge  orders. 
However,  estimates  were  made  of  the  approximate  order  number  and  the  results 
proved  to  be  satisfactory.  The  resultant  forces  between  the  disks  were  evaluated  by 
the  developed  program,  and  these  results  are  shown  in  the  Table  5.1.  Comparisons  of 
the  computed  contact  load  transfer  are  made  with  the  experimentally  measured 
values,  and  the  relative  errors  for  all  cases  are  less  than  10%. 
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5.4  Summary 

According  to  the  theoretical  results  and  experimental  data,  there  are  local 
differences  between  the  contact  stress  fields  of  cemented  disks  and  those  under 
Hertzian  dry-contact.  In  the  neighborhood  near  the  contacting  surface,  the  Hertzian 
contact  does  not  accurately  predict  the  proper  internal  stress  distribution.  Because  of 
this  a  new  hybrid  scheme  has  been  developed  using  non-Hertzian  contact  mechanics. 
This  technique  uses  a  family  of  contact  loading  polynomials  to  approximate  the  actual 
load  transfer  between  cemented  particles.  It  has  been  demonstrated  that  these  loading 
polynomials  produce  reasonably  accurate  stress  fields  which  match  with  experimental 
photoelastic  data,  and  can  be  incorporated  into  an  inverse  hybrid  scheme  whereby 
photoelastic  fringe  data  can  be  used  as  input  to  a  computer  code  in  order  to  determine 
the  interparticle  load  transfer.  The  polynomial  family  of  contact  stress  distributions 
can  cover  a  broad  variety  of  cement  types  of  various  moduli  ranging  from  soft  to 
hard.  The  accuracy  of  the  calculations  appeared  to  be  satisfactory  as  differences 
between  numerical  and  experimental  load  transfers  were  all  less  than  10%.  Further 
verification  of  this  technique  will  be  done  on  additional  data. 
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TABLE  5. 1  Comparison  of  the  Experimental  and  Numerical  Results 


Case  1  Urethane 

Case  4  Aluminum 

fSSt 

Num* 

Err. 

m 

Num* 

Err. 

Num* 

Err. 

Num* 

Err 

200 

192 

4.0% 

267 

268 

0.3% 

95 

100 

4.3% 

91 

97 

6.0% 

222 

in 

2.4% 

690 

642 

7.0% 

191 

204 

6.5% 

134 

138 

3.0% 

1291 

1293 

0.2% 

270 

293 

8.5% 

181 

191 

5.1% 

1389 

1387 

0.1% 

357 

380 

6.6% 

225 

221 

1.7% 

■ 

mi 

1811 

2.2% 

447 

481 

7.5% 

270 

280 

3.6% 

1905 

1968 

3.3% 

537 

573 

6.6% 

314 

332 

6.0% 

■ 

2226 

2367 

6.3% 

629 

646 

2.8% 

359 

378 

5.3% 

■ 

111 

751 

4.1% 

407 

422 

3.8% 

812 

855 

5.4% 

448 

453 

1.2% 

*  All  load  units  are  in  Newtons. 
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(c)  No  Cement  (Hertzian  Case) 

Fig  5  1  Photoelastic  Fringe  Patterns  of  Disks  in  Contact 
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(a)  Coordinate  System  and  Loading  Geometry 


(F=312  lb.) 


y/R 


(b)  Maximum  Shear  Stresses  Along  y-axis  of  the  Disk  Cemented  by  Super-Glue 
Fig.  5.2  Comparison  of  the  Maximum  Shear  Stresses 
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CHAPTER  6 


PRELIMINARY  STUDIES  OF  DYNAMIC  LOAD  TRANSFER 
THROUGH  SATURATED  GRANULAR  MEDIA 
6.1  Introduction 

Preliininary  microstructural  modeling  efforts  have  been  conducted  to  investigate  the 
effect  of  pore  fluid  in  granular  and  particulate  materials.  The  emphasis  in  this  segme  mr 
program  is  to  establish  reasonable  models  of  the  local  microstructural  interparticle  ^  lact 
effects  associated  with  pore  fluid,  and  to  relate  these  effects  to  the  wave  propagational 
characteristics  through  model  granular  materials. 

Our  efforts  in  modeling  the  dynamic  response  of  such  fluid-saturated  particulate 
materials  have  used  the  discrete  element  method.  This  computational  method  has  shown  to 
be  very  useful  to  model  the  response  of  discontinuous  and  particulate  materials  by  studying 
the  behavior  of  individual  idealized  particles  in  assembly  systems.  For  applications  to  wave 
propagation,  this  scheme  employs  large  assemblies  of  idealized  particles  to  model  the  dynamic 
behavior  of  granular  materials.  Numerous  studies  by  Sadd,  et.al.  (1989,  1991  and  1993)  have 
shown  that  such  modeling  techniques  can  simulate  granular  materials  and  predict  results 
which  compare  with  experimental  data.  The  numerical  strategy  uses  Newtonian  rigid-body 
dynamics  to  calculate  the  translational  and  rotational  motion  of  each  particle  in  these  model 
assemblies.  In  this  fashion,  the  dynamic/acoustic  response  of  the  model  system  may  be 
determined,  and  parameters  such  as  wave  speed  and  amplitude  attenuation  can  be  calculated 
for  specific  model  assemblies.  Furthermore,  these  wave  propagat-  al  characteristics  can  be 
related  to  the  model  material's  microstmcture  or  fabric. 


Of  fiindamoital  importance  for  the  proper  application  of  this  numerical  scheme  is  the 
development  of  correct  interparticle  contact  laws  which  predict  the  contact  forces  between 
adjacent  particles  in  the  assembly.  These  laws  play  a  key  role  in  determining  the  constitutive 
response  at  the  micro-level.  Much  of  our  beginning  work  focussed  on  developing  new 
contact  models  for  particles  interacting  with  each  other  through  a  fluid.  The  presence  of  pore 
fluid  can  have  significant  effects  on  the  dynamic  response  of  granular  materials  by  adding 
viscous  drag  forces  on  particles  and  by  changing  the  contact  response  between  adjacent 
particles  through  squeeze-film  elastohydrocfynamic  action. 

6.2  The  Discrete  Element  Method 

As  mentioned,  the  discrete  or  distinct  element  method  is  a  modeling  strategy  which 
uses  Newtonian  rigid-body  mechanics  to  model  the  translational  and  rotational  motion  of 
particles  in  model  material  assemblies  (Figure  6.1).  Contact  laws  between  adjacent  particles 
are  constructed  which  serve  to  determine  the  contact  force  as  a  function  of  the  relative 
displacement  or  relative  velocity  between  the  particles.  Applying  Newton's  law  to  the  i-th 
particle  would  yield 

n 

^  ^  (6.1) 

i-  I 

where  are  the  j-contact  forces  on  the  i-th  particle  and  represents  any  non-contacting 
forces  (e  g.  viscous  drag).  Equation  (6.1)  would  then  yield  the  particle  acceleration  with 
given  contact  and  non-contact  forces.  The  technique  establishes  a  discretized  time  stepping 
numerical  routine,  in  which  granule  velocities  and  positions  are  obtained  from  numerical 
integration  of  the  computed  accelerations.  It  is  assumed  that  during  each  time  step, 
disturbances  cannot  propagate  from  any  particle  further  than  its  immediate  neighbors.  Under 
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these  assumptions,  the  method  becomes  explicit,  and  therefore  at  any  time  increment  the 
resultant  forces  on  (and  thus  the  accelerations  of)  any  particle  are  determined  solely  by  its 
immediate  neighbor  interactions.  In  wave  propagation  applications,  the  movements  of  the 
individual  particles  are  a  result  of  the  propagation  through  the  medium  of  disturbances 
originating  at  particular  input  loading  points.  Consequently,  wave  speed  and  amplitude 
attenuation  (intergranular  contact  force)  will  be  functions  of  the  physical  properties  of  the 
discrete  medium,  i.e.  the  microstructure.  Figure  6.2  illustrates  a  typical  flowchart  of  the 
discrete  element  method.  In  order  to  model  particulate  materials,  the  method  is  typically 
applied  to  granular  systems  containing  large  numbers  of  idealized  particles  (e  g.  circular  disks 
or  spheres)  in  regular  or  random  packing  geometries.  Such  model  materials  can  be 
computationally  generated  with  varying  degrees  of  microstructural  fabric  and  pore  fluid 
content. 

6.3  Contact  Laws  for  Particles  in  Saturated  Media 

Our  modeling  efforts  began  with  the  development  of  new  contact  laws  governing 
particle  interactions  through  a  fluid.  As  mentioned,  the  presence  of  pore-fluid  will  create 
several  new  forces  at  the  particulate  level.  For  saturated  conditions,  capillary  forces  from  the 
pore  fluid  would  be  negligible;  whereas  the  phenomena  of  viscous  drag  and  squeeze-film 
contact  lubrication  could  produce  sizeable  particulate  forces.  Viscous  drag  effects  would  be 
predominant  for  large  particle  motions  with  significant  inter-particle  spacings  where 
neighboring  interaction  forces  are  vanishingly  small.  On  the  other  hand,  squeeze-film  contact 
forces  would  be  the  most  significant  for  dense  packings  of  particulate  systems,  and  this  case 
is  the  situation  to  be  studied. 

In  order  to  model  this  squeeze-film  effect,  consider  the  case  of  two  particles 
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(q)proxiinated  as  two  circular  disks)  onbedded  in  a  viscous  fluid,  and  approaching  each  other 
with  a  relative  normal  velocity  yfij,  see  Figure  6.3.  The  fluid  film  in  the  contact  zone  will  be 
squeezed,  and  a  sizeable  pressure  is  thus  built  up  in  this  contact  region.  The  fluid  pressure 
distribution  in  the  gap  will  exert  a  loading  on  each  particle  surface,  and  this  loading  is 
sufficient  in  certain  cases  to  produce  significant  particle  deformation.  In  this  fashion,  the 
contact  law  governing  particle  interactions  for  saturated  granular  materials  changes 
considerably  fi'om  the  contact  mechanics  for  dry  materials.  The  saturated  case  is  more 
complicated  since  it  involves  a  coupling  of  both  fluid  and  solid  behaviors. 

The  contact  force  for  such  a  problem  can  be  modeled  using  the  theory  of 
elastohycb-ocfynamics  e.g.,  Cameron  (1967)  and  Gohar  (1988).  The  schematic  of  the  two-disk 
squeeze-film  problem  shown  in  Figure  6.3  describes  the  problem  geometry.  The  distance 
between  particle  surfeces  is  denoted  by  h(x,t),  a  fimction  of  coordinate  x  and  time  t,  and  S(x,t) 
is  the  granule  deformation.  Under  certain  assumptions  the  fluid  pressure  p(x,()  can  be  related 
to  h(x,tj  by  the  Reynolds  equation  of  lubrication  theory 


±0,^lZ)  .  12„^ 

dx  dx  dt 


(6.2) 


where  t]  is  the  viscosity  of  the  fluid. 

Since  h(x,t)  can  include  particle  deformation  d(x,t)  which  is  determined  by  elasticity 
theory,  equation  (6.2)  is  a  coupled  relation,  and  typically  an  iterative  numerical  method  would 
be  required  to  solve  this  equation  for  the  pressure  distribution p(x,t).  Incorporating  such  a 
numerical  scheme  into  the  discrete  element  method  would  result  in  a  very  computationally 
intensive  procedure,  and  would  therefore  result  in  unreasonable  amounts  of  CPU  time  to  run 
the  simulations.  Therefore  a  simplified  alternative  procedure  was  developed  in  order  to  avoid 


the  numerical  complexities.  This  simplified  approach  first  assumes  that  the  particle  is  rigid, 
and  this  allows  a  closed-form  solution  for  p(x,i)  from  equation  (6.2).  Once  this  pressure 
distribution  is  found,  the  particle  deformation  may  be  calculated  from  simple  Hertz  contact 
stress  theory,  and  this  new  particle  shape  may  be  then  used  to  calculate  a  new  pressure 
distribution.  Thus  initially  let  S(x,t)  =  0,  and  the  undeformed  particle  surface  can  be 
approximated  by  the  relation 

X  ^ 

A(x./)  «  /<0.0  +  —  (6.3) 
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where  R  is  the  particle  radius.  The  solution  to  the  Reynolds  equation  (6.2)  for  this  case  is 
found  to  be 

p(xy)  =- 6Rii-^5L  (5  4) 

hHx/) 

With  the  pressure  distribution  known,  the  total  contact  load  between  the  particles  is  given  by 
integrating  the  pressure,  i.e.  F  =  H  J*  pdx,  where  H  is  the  disk  thickness.  Carrying  out  the 

integration  over  the  disk  surface  yields 

r  /  « 

m  =  3V2itTi/f  m  (6.5) 

I  mol 

In  order  to  include  the  effects  of  particle  deformation,  a  combined  model  with  two, 
series-connected  stiffiiess  (spring)  elements  may  be  used.  One  spring  element  represents  the 
fluid  stiffiiess  wWle  the  other  element  includes  the  solid  (particle)  stiffiiess  resulting  from  the 
elastic  deformation  within  the  particle.  The  pressure  loading  will  produce  deformations  such 
that  the  distance  between  particle  surfaces  will  change  from  h  to  h(x,t)  +  6(x,t),  where  6(x,t) 
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is  the  contact  deformation  for  the  dry  case.  If  the  disk  shi^)e  is  assumed  to  be  unchanged,  then 
the  new  load  F(i)  is  still  given  by  equation  (6.5)  using  the  modified  value  of  h(x,t),  i.e. 


where  C,  =  3'/2t]  tiR^.  Equation  (6.8)  then  gives  the  normal  contact  force  when  two  disks 
approach  each  other.  For  the  case  when  the  two  particles  move  apart,  contact  law  (6.8)  is 
not  suitable  since  a  negative  velocity  V,  may  lead  to  a  negative  value  within  the  square  root. 
For  this  case  the  solid  contact  law  will  be  us>  i 


The  above  results  address  only  the  normal  contact  response  between  particles.  For 
the  tangential  contact  behavior,  a  simple  modeling  scheme  based  upon  the  simple  shearing  of 


a  Newtonian  fluid  is  applied  within  the  small  interparticle  gap  as  shown  in  Figure  6.4.  If  the 
gap  is  smaU  then  the  fluid  motion  may  be  approximated  by  a  linearly  varying  distribution  (see 
Fig.  6.4a),  and  the  resulting  shear  stress  between  the  particles  will  then  be  given  by 


where  U  is  the  relative  tangential  velocity  between  particles.  Thus  the  force  between  any  two 
generic  particles  can  be  estimated  by  the  product  of  this  shear  stress  times  an  effective  area 
within  the  gap  zone.  We  are  currently  experimenting  with  various  schemes  to  determine 
appropriate  effective  areas  for  the  calculation  of  the  tangential  force. 

Although  this  simple  scheme  provides  reasonable  estimates,  it  will  break  down  if  the 
interparticle  gap  becomes  small.  For  this  case,  the  spacing  parameter  h  will  be  very  small,  and 
thus  the  velocity  gradient  and  hence  the  predicted  shear  stress  will  become  quite  large. 
Therefore  this  model  will  predict  unreasonably  high  values  of  the  interparticle  tangential  load 
transfer.  This  situation  has  been  studied  previously  by  the  lubrication  community,  and  it  has 
been  found  that  for  very  small  gap  widths,  the  assumption  of  zero  slip  conditions  is  not 
satisfied.  Previous  studies  have  found  that  it  is  necessary  to  introduce  a  slip  velocity,  to 
account  for  the  fact  that  fluid  particles  will  not  adhere  to  the  solid  particle  surfaces.  This 
situation  is  shown  in  Fig  6.4b,  and  the  effect  of  this  slip  velocity  is  to  lessen  the  severe 
velocity  gradient  within  the  gap.  This  modeling  approach  will  control  the  predicted  tangential 
load  transfer,  and  we  intend  to  formulate  this  scheme  into  the  wave  propagation  DEM  code 
6.4  Results 

In  order  to  investigate  the  merit  of  the  proposed  modeling  scheme,  several  one-  and 
two-dimensional  material  models  were  created,  and  the  discrete  element  computer  code  was 
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thoi  applied  to  each  of  these  models  under  specific  dynamic  loading  conditions  designed  to 
simulate  geoacoustic  wave  propagation. 

6. 4. 1  One-Dimensional,  Single  Chain  Simulations 

The  previously  described  contact  law  has  now  been  incorporated  in  our  basic  discrete 
elanent  computer  code.  Preliminary  computer  runs  have  been  conducted  in  order  to  simulate 
the  acoustic  wave  propagation  in  an  one-dimensional  granular  material  simulated  by  a  single 
strrught  chain  of  circular  disks  as  shown  in  Figure  6.5.  Considering  the  modeling  of  a  coarse 
sand,  the  circular  particles  had  a  diameter  of  2mm,  thickness  of  1mm,  and  were  assumed  to 
be  of  quartz  with  a  density  of  2650  kg/m^.  The  simulations  incorporated  different  gap 
openings  (of  the  order  of  one  miaometer)  between  each  neighboring  particles,  and  the  model 
system  was  assumed  to  be  completely  saturated  with  water  of  viscosity  of  0.0008  N-s/m^. 
The  input  loading  applied  to  one  end  of  the  particle  chain  was  taken  to  be  a  triangular  time 
dependent  pulse  of  duration  2.5ps  with  a  peak  value  of  4N  to  simulate  a  wave  front  acoustic 
pressure  of  10*  N/m^  (i.e.  10  bar). 

Using  these  parameters,  a  computer  simulation  was  run  and  the  wave  attenuation 
results  are  shown  in  Figure  6.6.  It  was  found  that  the  attenuation  is  quite  large  during  the 
passage  through  the  first  several  particles,  but  then  the  attenuation  rate  decreases  with  the 
propagational  distance.  These  preliminary  simulations  indicated  that  the  wave  speed  was 
approximately  840  m/s.  The  effects  of  the  gap  distance  and  viscosity  on  the  wave  speed  was 
also  investigated,  and  these  results  are  shown  in  Figure  6.7.  As  expected,  the  wave  speed 
decreases  with  the  gap  distance,  and  when  the  gap  is  small,  tlii  effects  are  rather 
pronounced.  An  increase  of  fluid  viscosity  also  leads  to  a  decrease  of  the  wave  speed,  and 
again  this  variation  was  most  severe  for  the  case  of  small  gap  distances 
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6.4.2  Two-Dimensional  Simulations 


In  order  to  investigate  the  use  of  the  model  for  a  multi-dimensional  simulation,  a  two- 
dimensional  granular  model  assembly  was  get^ated.  The  assembly,  shown  in  Figure  6.8  was 
created  using  one  of  our  computational  random  media  generators.  These  generating  codes 
can  create  large  random  particulate  assemblies  with  varying  degrees  of  microstructural  fabric 
such  as  porosity,  coordination  number,  particle  contact  normal  and  branch  vector 
distributions,  void  vectors,  etc.  Dynamic  input  for  wave  propagation  was  created  through 
simultaneous  loading  of  particles  along  left  and  bottom  edges  of  the  assembly  as  shown. 
Input  loadings  were  applied  with  identical  magnitude  and  time  history  as  used  in  the  one¬ 
dimensional  simulations.  The  transmitted  wave  output  (measured  by  the  inter-particle  contact 
forces)  was  collected  among  the  particles  along  the  right  and  top  sides  of  the  assembly,  and 
thus  this  technique  can  provide  comparisons  of  horizontal  and  vertical  wave  propagationaJ 
behaviors. 

In  order  to  correlate  material  microstructure  to  the  wave  propagational  behaviors,  the 
branch  vector  distribution  of  the  assembly  in  Fig.  6.8  was  computed.  The  branch  vector  is 
a  commonly  used  microstructural  fabric  measure  that  relates  the  relative  positions  between 
particles.  It  is  defined  as  the  vector  drawn  between  adjacent  particle  mass  centers,  and  Fig 
6  9  illustrates  a  polar  rose  diagram  plot  of  this  distribution  for  the  assembly  of  Fig  6.8.  It  can 
be  observed  from  the  distribution  plot  that  this  model  material  system  shows  preferred 
directions  (anisotropy)  with  more  branch  vectors  distributed  in  the  horizontal  direction  than 
in  !hr  vertical.  Thus  we  would  expect  more  wave  transmission  for  the  horizontal  loading 
case 

Discrete  element  simulation  results  of  the  wave  propagation  through  this  assembly  are 
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shown  in  Figure  6. 10.  Normalized  average  wave  transmission  through  the  assembly  for  the 
horizontal  and  vertical  loading  cases  are  illustrated.  The  vertical  scale  represents  the  average 
of  all  interparticle  load  transfers  at  the  output  s*de  of  the  assembly  normalized  with  respect 
to  the  input.  These  results  indicate  that  the  wave  transmission  is  clearly  related  to  the 
anisotropic  fabric  of  the  model  particulate  material  in  that  less  wave  attenuation  occurs  for 
the  horizontal  propagation  case.  Figures  6.11  and  6.12  show  vertical  and  horizontal 
propagation  results  of  the  same  model  material  for  the  case  with  different  inter-particle  gap 
spacing.  It  can  be  observed  that  with  an  increase  of  the  gap  spacing,  attenuation  will  increase 
while  propagational  speed  will  decrease.  These  results  are  consistent  with  our  one¬ 
dimensional  simulations. 

Future  work  with  saturated  systems  will  include  improvements  in  the  basic  inter- 
particle  modeling  and  we  mtend  to  conduct  an  e^  tensive  series  of  computer  runs  to  establish 
more  specific  relationships  with  fabric  and  wave  propagational  behaviors. 
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Fig.  6.3  Elastohydrodynamic  Action  Between  Particles 
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Fig.  6.5  One-Dimensional  Single  Chain  Model 
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6.6  Wave  Attenualion  Behavior  for  Sinule  Chain  Model 
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Fig.  6.7  Wave  Speed  as  a  Function  of  Gap  Spacing 
and  Viscosity  for  Single  Chain  Model 


Fig  6  8  I'wo- Dimensional  Random  (iramilar  Model 


Fig.  6.9  Branch  Vector  Distribution  Diagram  of  Random  Model 
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CHAPTER  7 

*  INFLUENCE  OF  PARTICLE  SHAPE  ON  THE 

DYNAMIC  LOAD  TRANSFER  PROCESS 

A  preliminary  investigation  was  conducted  to  study  the  effect  of  contact  proHle 
and  granular  media  porosity  on  inter-granular  load  transfer,  wave  velocity,  and 
^  wavelength  of  a  stress  wave  propagating  through  the  granular  media  subjected  to 

explosive  loading.  Dynamic  photoelasticity,  the  combination  of  high  speed  photography 
and  photoelasticity,  was  used  to  study  stress  wave  transmission  through  the  granular 
§  media  that  were  simulated  by  assemblies  of  circular  or  non-circular  particles  in  contact. 

TTie  photographs  thus  obtained  were  analyzed  to  determine  contact  load  as  a  function  of 
time,  wave  velocity,  and  wavelength.  Specimens  were  made  of  Homalite-100,  an  optically 
§  birefringent  material.  Since  some  of  the  specimens  were  not  circularly  symmetric, 

experiments  were  conducted  with  the  specimens  oriented  both  with  the  major  physical 
feature  toward  the  incoming  pulse  and  away  from  the  incoming  pulse.  The  results  indicate 
®  that  contact  profile  and  media  porosity  have  an  effect  on  wave  velocity  and  wavelength. 

7.1  Introduction 

^  Several  investigations  have  previously  been  conducted  by  many  researchers  to 

study  the  dynamic  behavior  of  granular  media.  The  load  transfer  characteristics,  wave 
velocity,  and  wave  dispersion  for  various  loading  parameters  and  microstructure  of 
I  granular  media  have  been  studied  in  detail.  The  effect  of  particle  shape  on  wave 

propagation  has  been  studied  by  Shukla  et.al  (1991)  using  elliptic  particles.  They  showed 
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that  the  wave  velocity  was  dependent  on  the  particle  shape.  The  wavelength  of  the  stress 
wave  pulse  showed  no  change  with  the  particle  shape  as  long  as  the  contact  interval  ( the 
number  of  contact  points  per  unit  wavelength  )  is  kept  constant.  The  load  transfer 
characteristics  change  with  variations  in  particle  shape.  However,  in  most  of  the  past  work 
the  granular  media  was  simulated  by  assemblies  of  circular  or  elliptical  discs.  All  contacts 
were  therefore  between  two  convex  surfaces.  In  real  life,  some  contacts  are  between 
concave  and  convex  surface.  The  change  of  particle  shape  and  thus  of  contact  profile 
causes  the  contact  stiffness  to  change.  If  the  radius  of  curvature  at  the  contact  is 
increased,  while  keeping  the  contact  interval  constant,  the  contact  stiffoess  will  increase, 
and  this  would  imply  a  higher  wave  velocity.  This  may  also  change  the  wave 
propagation  characteristics. 

In  this  work,  four  groups  of  disks  were  used  to  simulate  the  granular  media.  The 
interest  lay  in  the  investigation  of  the  effect  of  a  concave-convex  contact  and  porosity  of 
media  on  the  load  transfer  characteristics,  wave  speed  and  wavelength.  Specimen 
geometric  parameters  are  shown  in  Table  7.1.  The  Groups  1  and  2,  as  shown  in  Fig.  7.1, 
were  also  used  to  study  the  effect  of  incident  direction  on  wave  propagation 
characteristics. 

7.2  Experimental  and  analytical  procedure 

The  method  of  dynamic  photoelasticity  was  used  as  the  experimental  technique 
throughout  this  investigation.  This  method  has  been  shown  by  many  researchers  to  be 
very  efficient  in  the  study  of  dynamic  wave  characteristics  in  granular  media  (for  example 
Shukla  1991).  As  stated  earlier,  Homalite-100  was  chosen  to  make  specimens.  To 
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fabricate  the  various  particles,  one  inch  diameter  circular  disks  were  machined  into  the 
geometries  shown  in  Table  7.1.  A  cooling  system  was  used  during  machining  of  the 
specimens  to  avoid  producing  any  heat  related  residual  stresses  in  the  specimens. 

The  granular  media  was  simulated  by  single  chain  assemblies  of  similarly  shaped 
particles  in  contact.  One  such  setup  is  shown  in  Fig.  7.1  (a).  The  same  shaped  particles 
were  also  used  in  an  alternative  setup,  as  shown  in  Fig.  7.1  (b),  to  investigate  the  effect 
of  incident  direction  of  the  wave  pulse  on  the  propagation  characteristics.  The  six 
different  single  chain  assemblies  used  in  this  study  are  shown  in  Fig.7.2. 

The  experimental  models  were  placed  in  the  optical  bench  of  a  high  speed  Cranz- 
Schardin  type  camera.  The  particles  were  loaded  by  exploding  a  small  amount  (10  mg) 
of  Lead  Azide  in  a  specially  designed  charge  holder.  The  camera  was  triggered  at  some 
prescribed  delay  time  after  igniting  the  explosive  that  was  placed  on  top  of  the  model  (see 
Fig.  7.1).  The  high  speed  photography  system  operated  as  a  series  of  high  intensity, 
extremely  short  duration,  pulses  of  light  and  provided  twenty  photoelastic  images  at 
discrete  time  during  the  dynamic  event. 

A  typical  sequence  of  six  images  for  the  single  chain  setup  is  shown  in  Fig.  7.7. 
These  images  of  the  propagation  phenomenon  were  enlarged  and  digitized  to  facilitate 
analysis.  The  Hertz  stress  field  equations  along  with  the  stress  optic  law  (see  Chapter  4) 
were  used  to  calculate  the  contact  load  at  a  particular  contact  point  and  time  by  applying 
the  multi-point,  non-linear  least  square  method  suggested  by  Shukla  and  Nigam(1986). 
The  wavelength  of  the  stress  wave  pulse  can  be  measured  directly  from  the  photographs 
of  the  wave  propagation  process  (Fig.  7.7).  The  propagation  distance  is  plotted  as  a 
function  of  time  to  obtain  the  wave  velocity. 
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7  J  Results  and  discussion 


To  understand  the  effect  of  particle  shape  on  wave  velocity  ,  wavelength  and  load 
transfer  characteristics,  a  series  of  experiments  were  conducted  with  the  experimental 
arrangement  shown  in  Fig.  7.1. 

Figs.  7.7  and  7.3  show  typical  isochromatic  fringe  patterns  and  average  velocity 
of  the  stress  pulse  for  chain  1  (see  Fig.  7.2)  respectively.  Likewise,  the  fringe  patterns 
and  average  velocity  for  chains  2,  3  and  4  are  shown  in  Figs.  7.8,  7.9,  7.10,  7.4,  7.5,  and 
7.6. 

The  average  velocities  of  stress  pulse  for  chains  1  and  2  ( contact  between  convex 
and  flat  surfaces  )  were  determined  to  be  approximately  equal,  1010  m/s  and  1020  m/s 
respectively  (  see  Figs.  7.3  and  7.4  ).  This  shows  that  the  incident  direction  of 
propagation  has  no  effect  on  average  pulse  velocity  for  this  particle  shape.  Figs.  7.5  and 
7.6  show  that  the  average  velocities  of  stress  pulse  for  chain  3  and  4  (  contact  between 
convex  and  concave  surfaces)  are  1140  m/s  and  1100  m/s  respectively.  This  indicates 
that  the  incident  direction  of  propagation  also  has  no  effect  on  average  pulse  velocity  for 
this  particle  shape.  The  stress  pulse  velocity  in  chain  3  and  4  is  10%  higher  than  that  in 
chain  1  and  2.  This  can  be  explained  by  the  fact  that  the  contact  stiffness  between  convex 
and  concave  surfaces  is  higher  than  that  between  convex  and  flat  surfaces. 

Compared  with  the  pulse  velocity  in  one  inch  diameter  circular  disk  chain,  the 
wave  velocity  in  chain  3  and  4  is  slightly  higher  because  of  larger  contact  stiffness. 
However,  the  wave  velocity  in  chain  1  and  chain  2  is  slightly  lower  than  that  in  one  inch 
diameter  circular  disk  chain  despite  the  larger  contact  stifftiess.  The  reason  for  this  is  the 
increased  number  of  contacts  per  unit  length  in  chains  1  and  2  as  compared  with  a  one 
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inch  diameter  circular  disk  chain. 


Typical  isochromatic  fringes  obtained  from  various  experiments  involving  single 
chain  assemblies  are  shown  in  Fig.s  7.7,  7.8,  7.9,  and  7.10.  The  wavelengths  for  chains 
1,2,3,  and  4  were  directly  measured  from  these  Figs,  and  are  shown  in  Table  7.2.  The 
pulse  lengths  for  chains  1  and  2  were  found  to  be  equal.  Chains  3  and  4  also  produced 
equal  pulse  lengths.  This  would  lead  one  to  believe  that  the  incident  direction  of 
propagation  has  no  effect  on  the  pulse  length.  It  was  also  noticed  that  wavelengths  of 
chains  1  and  2  were  larger  than  those  of  chains  3  and  4.  The  pulse  lengths  of  chains  3 
and  4  were  nearly  four  disk  lengths.  This  difference  in  pulse  lengths  may  be  the  effect 
of  contact  profile.  Fig.  7.10  also  shows  that  the  fringe  patterns  tend  to  take  the  shape  of 
the  particle  as  the  stress  wave  travels  through  it.  This  phenomenon  also  can  be  seen  in 
elliptical  disk  experiments.  (Shukla  1991) 

Data  collected  from  fringe  patterns  in  Figs.  7.7,  7.8,  7.9,  and  7.10.  was  used  to 
calculate  contact  load.  The  contact  load  variations  with  respect  to  time,  for  different  disk 
chains  are  shown  in  Figs.  7.11,  7.12,  7.13,  and  7.14.  The  contact  loads  have  been 
normalized  with  respect  to  the  peak  contact  load  between  particles  5  and  6.  This 
normalization  allows  for  easy  comparison  between  different  experiments. 

Normalized  contact  load  plots  provide  information  about  load  attenuation  as  the 
stress  wave  travels  from  contact  to  contact.  They  also  provide  information  about  the 
duration  of  contact.  Figs.  7.11,  7.12  and  Figs.  7.13,  7.14  show  that  load  duration  for  flat 
particles  was  slightly  longer  than  that  for  concave  particles.  The  reast^n  is  that  the  concave 
particle  chain  has  higher  wave  velocity  and  shorter  wavelength  than  flat  particle  chain 
has. 
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Figs.  7.15  and  7.16  are  the  average  velocity  plots  for  chains  5  and  6.  If  we  define 
the  porosity,  t),  as  in  equation  (7.1), 

area  of  hole  „  i\ 

area  of  disk 

the  velocity  plots  show  that  when  r\  increases,  i.e.  the  hole  is  made  larger,  the  average 
velocity  decreases.  Compared  with  the  wave  velocity  in  circular  disks,  the  average 
velocity  changes  very  little  for  small  values  of  q,  i.e.  a  very  small  hole.  It  appears  that 
there  may  be  a  threshold  at  which  the  size  of  the  hole  becomes  a  significant  factor  in  the 
pulse  velocity. 

From  Figs.  7.17  and  7.18,  the  fringe  pattern  becomes  distorted  because  of  the 
boundary  effects  produced  by  the  hole.  When  the  pulse  was  first  entering  a  disk  with  a 
hole,  and  the  contact  load  was  still  relatively  small,  the  fringes  appear  Hertzian  in 
shape(see  Fig.  7.17,  Frame  9).  As  the  load  increases,  the  fringe  patterns  show  the 
boundary  effects(see  Fig.  7.16,  Frame  13).  Because  of  the  boundary  effects.  Hertz  contact 
equations  cannot  be  used  to  analyze  the  data. 

Fringe  patterns  in  Figs.  7.17  and  7.18  show  that  a  lot  of  energy  remains  in  the 
disks  after  the  main  pulse  has  already  passed  through.  Photoelastic  experiments  cannot 
be  used  to  measure  the  wavelength  exactly  due  to  the  limited  field  of  view  provided  by 
the  optical  equipment. 

7.4  Conclusions 

The  results  indicate  that: 
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•  The  contact  profile  has  some  effect  on  wave  velocity  and  wavelength.  But. 
the  difference  between  the  wave  velocities  for  various  geometries  was  only 
10%.  Therefore,  more  data  are  needed  to  make  a  conclusive  statement. 

•  For  the  same  contact  profile,  the  incident  direction  of  propagation  does  not 
appear  to  change  the  wavelength. 

•  The  porosity  of  granular  media  changes  the  wave  velocity.  The  wave 
velocity  decreases  as  porosity  increase. 

Most  importantly,  the  experiments  that  have  been  conducted  in  this  area  are  in  the 
preliminary  stages  and  further  study  needs  to  be  done  to  obtain  any  conclusive  results. 
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Table  7.1  Geometric  parameters  for  different  particles  used  in  the  experiments. 


Group  No.  Specimen  Rl/Dl  (inch)  R2/D2  (inch)  h  (inch) 


Geometry 


Table  7.2  Experimental  result  of  wave  velocity  and  wavelength  for  concave  and  flat 
particle 
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Fig.  7, 2  The  six  different  single  chain  assemblies  j  ■ 
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Rg.7.3  Wave  propagation  distance  versus  time  for  chain  1  in  fig.7.2 
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Fig. 7. 4  Wave  propagation  distance  versus  time  for  cham  2  in  fig.7.2 
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Fig,7.5  Wave  propagation  distance  versus  time  for  chain  3  in  fig.7.2 
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Fig.7.6  Wave  propagation  distance  versus  time  for  chain  4  in  fig.7.2 
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Fig.7.8  Typical  isochromatic  fringes  obtained  for  a  single  chain  experiment  eondue 
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Fig.7.10  Typical  isochromatic  fringes  obtained  for  a  single  chain  experiment  conducted 

using  chain  4  fig. 7.2 
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Fig.7.13  Normalized  contact  load  as  a  function  of  time  for  chain  3  in  fig.7.2 
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Fig. 7.14  Normalized  contact  load  as  a  function  of  time  for  chain  4  in  fig.7.2 
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Fig.7.15  Wave  propagation  distance  versus  time  for  chain  5  in  fig.7.2 
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Fig.7.16  Wave  propagation  distance  versus  time  for  cha'n  6  in  fig.7.2 
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